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FOREWORD 


This report, prepared by the Dynamics and Loads' Section, Martin 
Marietta Corporation, Denver Division, under Contract NAS8-29944, presents 
the technical approach and the results of a study contract for the dynamic 
characteristics of a spherical tank/fluid/bladder to be used in Solar 
Electric Propulsion Stage (SEPS) . The study was administered by the 
National Aeronautics and Space Administrations, George C. Marshall Space 
Flight Center, Huntsville, Alabama, under the direction of Mr. Frank Bugg, 
Systems Dynamics Laboratory. 
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INTRODUCTION 


The proposed Solar Electric Propulsion Stage (SEPS) will be 
responsible for the transmittal of Shuttle-based payloads between low 
orbit and geosynchronous orbit. The SEPS will be a reuseable vehicle 
powered by accelerated mercury ions. The mercury propellant will be 
contained in spherical tanks. Propellant orientation within the tank 
will be controlled by a neoprene expulsion bladder; propellant feed to 
the ion engines will be maintained by insertion of freon pressurant 
into the ullage space between the tank and the expulsion bladder. The 
mass of the mercury will be a very significant portion of the total mass 
of the spacecraft. It is, therefore, apparent that this high mass frac- 
tion cannot be ignored and its dynamic characteristics must be investigated 
This study obtained the equilibrium shapes, vibration modal properties 
and mechanical equivalent slosh models for the mercury and bladder for five 
fill conditions and two acceleration levels. 
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1 . SCOPE 

1 . 1 Purpose 

The purpose of this report is to document the investigation performed 
under contract NAS8-29944, "Finite Element Analysis of Mercury Slosh in 
the Solar Electric Propulsion Stage". 

1 . 2 Sco pe 

This report documents the dynamic characteristics of the system of 
spherical tank/hemispherical expulsion bladder/mercury propellant. The 
static equilibrium shape corresponding to various ullage and gravity con- 
figurations was established. The finite— element approach was ap 
plied to different fill conditions in different gravity fields to evaluate 
the lateral sloshing mode shapes and frequencies. The resulting mode 
shapes and frequencies was used to define a spring-mass mechanical 
analog that describes the sloshing phenomenon. Computer programs for the 
equilibrium shape and vibration analysis including the modeling were 
developed. 

1.3 Summar y 

The static equilibrium shapes of the neoprene bladder have been es- 
tablished corresponding to various ullage and gravity configurations under 
specified boundary conditions. The hemispherical bladder is taken to be 
attached at the diametral plane of the sphere with zero relative slope. 
With these shapes, the spherical tank with bladder and mercury has been 
modeled as an assemblage of finite-elements. The properties of these ele- 
ments have then been calculated using a linear displacement field. The 
dynamic charac teristics were obtained to be used to define a mechanical 
analog which will reproduce the sloshing phenomenon of the system. The 
computer programs for the static free surface and vibration analysis have 
been checked out. 

2. TECHNICAL APPROACH 

The problem was approached in two distinct steps. 

(1) Establish the static equilibrium shape. 

(2) Using the free surface shape from (1), model the system as an 
assemblage of finite elements and obtain the mode shapes and 
frequencies. Further define a mechanical analog for the system. 
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2.1 Static Equilibrium Shape - The method comprises of deter- 
mining the total energy of the system in terms of the variables involved 
and minimizing the total energy for the stable equilibrium state under 
the given constraints. This is outlined below. However, it is to be 
pointed out here that the idealized case here is for an infinitely long 
cylinder and the feeling being that the static free surface for this 
case will be close to the actual case of spherical tank. 


* 
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Introducing w(y) in the expression 




x(y)= >< = ±L R "> a J +R '‘ 4 * r ^'*) 


(4) 


The plus sign represents the upper point of the sphere and the negative 
sign the bottom and hence , 


"ft) ® 


Va 


(5) 


An assumption is made here that w(y) may be expressed in terms of a 
power series : 


M £ 

Substituting for w(y) from equation (6) in equation (5), 


(6) 


(7) 


The stretching in the membrane strip of length s measured along the 
curved shape is iT(s) . It is assumed that 


u<£>® 6* 


( 8 ) 


Therefore , 


<?> u-, —- . 6 

f'bS 



That is, B is the strain in the membrane. 

The total potential of the assumed configuration will now be calcu- 
lated. They are: 

(i) Bladder strain energy in bending: 




P* 


( 10 ) 


where: dV^ - strain energy in bending 

E = modulus of elasticity for the bladder 



X ■ moment of inertia 
dS = length 

p - radius of curvature 

Further 

p iw^rr 


Therefore, equation (10) can be rewritten as 


~ J_ El 
li p 

It can easily be seen that 

ifM 2 _ 3 a, 

nnn 

(12) 

r 1 

i 

i 

i 

(13) 

l 

f 

\ 

i 

Y- 

7N 

< 

ev- 

Substituting from (13) in (12) 



av = ±ei 
ft z 

( dfe Y' 

(14) 

/z 1 

Therefore, 



\= 
ft 2.4 

j 

o 

r \ ^ 

[-(%T] 6,i 

(15) 


where: h = thickness of the bladder 

(ii) The strain energy dud to membrane action: 


Resulting in 


dVj^ - ^(stress)^strain) h dS 


V - 

1*1 


[ (-, + (^f] ^ 


( 16 ) 


where V M = strain energy due to membrane action 
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(iii) Fluid gravitational potential energy: 

P H ’ye*-- *■*■)/* i 

Simplifying, 

V V 

v V z Pi f xV j ’ i 

6 O 

where: mass density of the fluid 

g = acceleration due to gravity 
xg 55 as shown in Figure 1 
V 0 = gravitational potential energy 

(iv) Work done by ullage pressure: 

dW = p.dV 


where: p “ ullage pressure 

Therefore, virtual work by ullage pressure due to virtual change 


o 

where: £>w = virtual change w 

- virtual work 

Substituting for &UT from (6) in (20), 


M 3°“ ; 


* 


as 


W = W (sq , a , a 2 5 ••**••**> ^ jij ) 


£ w - y r ^ w ^ w £6 


•i 


(17) 

(18) 


(19) 

in w, 

( 20 ) 

( 21 ) 

( 22 ) 

(23) 
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The comparison of equations (21) and (23) reveals. 


Cdc<i 


-H ‘ 


3 


for i=0, 1, 2, N 


and' 


(24) 


- O 

Cb B 


(25) 


In the present case, the system may be treated as conservative for a 
particular time, and therefore 

V = -W (26) 


where: V - is the total potential energy of the system 

W - work function 
Substituting in equation (24), 


and 


f 7>v 



^ _ 0 

r2>6 


(27a) 


(27b) 


under the boundary conditions of 

c 

(28) 

%S°) =° 

where : 2K = given ullage volume 

and 


V - V B + V M + V g 


(29) 


Assembling the expression 

ds -+ 


V- 

ZA 


/A t 


0<%T] s/z 


<i) & 


vz. 


i-?f c>-<) 




>(30) 
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The equation (30) has to be expressed in terms of a^'s and 

' X (>} - ^Cl s ') -- - Z <U f 


“ = ■ > ( R -f ) y - is * i 


XX d-XL 
XL =. — 


- - (** i ) - { ( **- il) - £ c(x ~ 0 > 


VF 

Substituting from (31), (32) and (33) in equation (30) and differentiating 
with respect to a^'s and b respectively and comparing with (27a) and 
(27b), 


3 *** 2 2 

?L&a* t 1 +W) 7 * 1 

JHa 




™L=0 = ^ 


for i = 0,1, 2 , 


V ZT* ^7*7 ^ 


Under the auxiliary Conditions of 

,J 7y*) * 

Ld'(o) - O 
j ter k 


ty (34a) 


(34b) 
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Finally, the equations of the system are 

£4* {-tfr-d'f't-m} frff Vi {-****/) 

u l IT , k 

+ c{i+&)}{-*j i-^*- 1 }] <*? 

+ n fwt-fm + = o 

and 

, .7* . 

** 6 /n+ {*jj ^-o 


for i = 0, 1, 2 , 

under the conditions 

1J (7* ) - o 

M'(°) ~ o 


N 


r 


it el <-t t fc 


7 


Thus , 


A= F fr,>Zf—>h+) 

In general, the equations can be written as, 


F 

'rv 


(%l)~ £Hh&£ . 

° Ct+{*?]K 

+ !%+(.*: ff-M i -M T*B*} 

c 

+ n f¥i*}{-f ri }J d s + t/fox * o 

o ▼ 


for n = 1, 2, 


, N+l 








(36) 


(37) 


(38) 


(39) 


V 

(40 a) 


(40b) 



10 


for n = N+2 


(40c) 


under the auxiliary conditions 

«r*o =° 

ca' (o) - o 
K 

The system of equations (40), i.e., (40a), (40b) are to be solved under 
conditions of (40c) for with iterative approach. Denoting (F n ^,, 
where s stands for iteration step, (F n ) stands for initial assumption, 




(41) 


In general, therefore, 


(% r 


(42) 


where n - 1, 2, N+2 and m = 1, 2, ..., N+2 (separately) 

Rewritting the equations (42) in matrix form: 

^Fi 


(F, 

Fz 

F3 


n/t, 




Ffjiz 


% 




rbFn+Z | ' btN+l 

Wi 1 'M* 






4- 


~4~ 


^3 



(43) 


The various terms involved in the matrix are as follows: 


(44a) 
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(c-^ac i+otff’V* otfotfl ^ 

C-{'-‘)OJ nflfc&iz-Ti) *? 

for n = 1, 2, ...... N+l and m = 1, 2, , N+l (separately) 


and 




-4 


(f f'^7 >7 ^ «*» 


for n = 1, 2, N+l and m ” N+2 (separately) 


Also 


{ 1I^: - £^£ 
n>%^ ' 




(44c) 


for n = N+2 and m = 1, 2, ...... N+l 

and 



for n = N+2 and m = N+2 

Recalling equation (13), the arc length can be written as 


(44d) 



(45) 


Assuming the bladder to be attached at the diametral plane, the original 
length of the arc is 

L ~ TT (46) 

Therefore, b can, now, be defined as 

6 = -^- (47) 

Thus, the system of equations (40a), (40b), (40c) and (47) are to be 
solved simultaneously with an iterative approach. This is accomplished 
in the computer program; 
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2.2 Vibration Analysis - Vibration analysis of structures having 
a large number of degrees of freedom is an ever present problem. Digital 
computer oriented techniques are primarily restricted by core size and 
computer time. Consequently, economically feasible eigenvalue/eigenvector 
techniques are needed for large size structural systems. 

However, for a large size structural system, if one is interested 
in a relatively small number of modes, Rayleigh- Ritz technique appears to 
have the advantage over the others. In this method, essentially a large 
problem is reduced to a smaller problem for the range of interest only. 

The solution is obtained as a linear combination of assumed linearly in- 
dependent mode shapes. 

The method consists in using the assumed mode shapes initially, 
which reduce the number of generalized coordinates used and then, by an 
iterative procedure, these modes are improved until they converge to the 
normal vibration modes of the structure. 

2.2.1 Rayleigh-Ritz Technique - The iterative Rayleigh-Ritz 
method (1) is used to calculate the mode shapes and frequencies of the 
system. This method is based on repeated application of the well known 
Rayleigh-Ritz technique using improved mode shapes for each iteration. 

This technique reduces the size of the system without degrading accuracy 
in the desired frequency range. The technique is briefly described 
here . 

For a discrete coordinate model of a structure having n degrees of 
freedom, the equations of motion can be written as 


["FHWH-W 


(4a) 


where : 


{*} = {*(*)} vector of discrete coordinate displacements, 
[M] = mass matrix 

= stiffness matrix 

If a solution of the type ^ j 


implying a simple 


harmonic motion is assumed, equations (48) can be written as 
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(CK3 - {*} = £°} (49) 

Equation (49) is recognized as a matrix eigenvalue problem of order n, 
whose eigenvectors M are the mode shapes and whose eigenvalues r* 2 j 
are the frequencies. A complete sequence of trial vectors 

{*} 3 > {*}*. ( 5 °) 

which are linearly independent, is assumed. The displacement {*} is 
then expressed as a linear sum of the first "m" trial vector, that is, 

H =[i]W <so 

(nx»n) (trtffj 

Substitution of equation (51) into (49) and for multiplying by 00^ gives 




(52) 


where : 


and 


C&1=KJ WM 
CH] =EfMKI 


(53) 

(54) 


Equation (52) is a matrix eigenvalue problem of reduced order "m" whose 
eigenvectors are cV) and eigenvalues are . The solution of 
equation (52) has the form 


{%} = crm*} 


(55) 


where : is the normalized coordinate vector. The eigenvalues, 

m , approximate the first "m" eigenvalues of the original struc- 
ture. The associated approximate eigenvectors DO of the original 
structure are obtained by substitution of equation (55) into (51), yielding 
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frwl) 

where : C*3 « C*3 tfl (57) 

The accuracy of the mode shapes Cv3 and frequencies Z&3 ob- 
tained depends entirely upon the trial vector in . If contains 

the true modal patterns, then the eigensolution for£*^^ and C«^3 are 
exact. However, in general, that is not the case. Exact results can 
be obtained for the first "m" modes of the structure if the trial vec- 
tors C*3 do not have any contribution from modes higher than "m". 

Thus, an improved set of trial vectors can be calculated by suppressing 
the contribution of higher modes in approximate mode shapes. The pro- 
cedure for suppressing the contribution of the higher modes is well 
known; in fact, it is the basis of the Power or Stodola-Vianello 1 matrix 
iteration method^ of modal analysis. Here, however, the method is 
applied to all modes simul taneously and is given as, 

[_<] [O - M M (58) 

The solution is carried out for [z] , which is then used to repeat 
equations (52) through (58) . The cycle can be repeated until all the 
mode shapes and frequencies CJ ^ have converged to within a pre- 
scribed tolerance. Convergence is assured because the technique is 
equivalent to a power iteration applied simultaneously to all modes. 

Thus, the convergence theorems associated with the power method are 
directly applicable. The role of the eigensolution (equation (52) 
is to prevent all modes from converging on the lowest mode. 

Associated with the iterative Rayleigh- Ritz technique are para- 
meters that affect the convergence and hence computer time which will 
be briefly discussed here. They are: 

(i) the initial mode shapes assumed to start the iteration process, 

(ii) the number of modes used, 

(iii) the repression of higher modes, and 

(iv) shifting. 
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(i) Initially Assumed Mode Sha p es - The choice of initial mode 
shapes plays a very important role in the success of the technique. 
Inherent with the initial mode shape selection are two basic problems * 

T 

(1) modes may be missed, and (2) the triple product QrqJ - [zj 
may be ill-conditioned if the columns of \Z] are not sufficiently inde- 
pendent. It does not appear that there is a way to guarantee that the 
above two conditions will be met with any selection of [ _l j , however, 
the chance of them occuring can be minimized with some judicious selec- 
tion of the vectors. Without proof or discussion, it is to be pointed 
out here that if the elements of the vector or of matrix [Z] are randomly 
generated, it has been found that the chances of the above two conditions 
being violated is very remote. 

(ii) Number of Modes Used - An increase in the number of modes 
used will, in general, decrease the number of iterations required for 
convergence. However, if more modes are used, the computer time for 
each iteration will increase because of the increase in sizes of the 
matrices used; hence, there is a tradeoff. 

(iii) Repression of Higher Modes - As pointed out earlier, exact 
results can be obtained for the first ,, m n inodes of the structure if the 
trial vectors in DO do not contain any contribution from modes higher 
than n m n . Generalizing, it can be said that an improved set of trial 
vectors can be calculated by suppressing the contribution from the higher 
modes in the approximate mode shapes at each step. This is achieved as 
follows . 

M. = [KfM iXL 

d * 

The subscript j denotes the iteration number. If this iteration is ' 

repeated sufficient number of times, modes corresponding to the lowest 

frequency will be reached. If this iteration is repeated too many times, 

the mode will repeat itself in one or more columns of [Z^ and will render 
T 

Z 3 L- 1 U] to be ill-conditioned. 

Its use here is not to converge to a mode but just to repress the 
higher modes and, hence, just a one time application is advisable. 
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However, in this case £kJ“ 1 is required for which £k] has to be 
non- singular. Thus, the technique can be applied only if £kJ is not 
singular or has been made such with some technique as described next. 

(iv) Shiftin g - Shifting is an useful technique to speed the con 1 
vergence of modes whose eigenvalues are close to the shift value. As 
an additional benefit of shifting process is the conversion of the 
stiffness matrix (in case of a free-free structure) from singular to a 
non- singular matrix. The method is as follows. 

The eigenvalue problem is 




(60) 


where: different quantities carry their usual meaning. 

Also, it has to be noticed that |^Kj may be singular. To introduce the 
shift value, , the following operation is performed. The quantity, 

A CM] i s su ^ tracte< ^ from both sides of equation (60) . Thus 


(m- a,w) ffl = M ffl 


(61) 


By definition 


and 


[SJ=[K> M") 


(62) 


(63) 


Therefore, final equation is 




(64) 


This is now . the eigen-problem to be solved rather than (60). It is to be 
noticed that is non- singular even if [k] was not. 

The eigenvalues of the original system are easily obtained as 

4=n*+A 5 

*** * 


( 65 ) 
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The convergence will be to the lowest absolute value of 2- . Thus, 
shifting by a value, , the eigenvalues, Gl> , around this shift 

point are converged to first. 

Some general remarks on Shift : 

(a) Analysis of a Free- Structure - Because a free structure has 

a singular stiffness matrix, the solution of the simultaneous 
equations in the iteration loop is not possible. However, the 
shift technique alleviates the problem. 

(b) Specific Frequency Range - When a shift value is used, the 
modes with eigenvalues closest to the shift value will con- 
verge first, which enables one to obtain the modes in the 
desired frequency range only. 

(c) Large number of modes - By repeated use of different shift 
values, any number of modes can be obtained. 

(d) The following observations are made without discussion. 

(i) If the lowest eigenvalues in the range , 6) 2 , , ec£, 

are needed, a shift value of zero should be used for a restrained 
structure and one for a free- free structure. 

(ii) If the modes are needed in an intermediate range, a shift 
midway between the lowest and the highest expected eigen- 
values should be used. 

2.2.2 Mass and Stiffness Matrices - The iterative Rayleigh- Ritz 
analysis subroutines require as input mass and stiffness matrices of a 
structure. To reduce engineering time required to perform an analysis, 
subroutines were included to calculate mass and stiffness matrices for 
general standard structural elements. The basic idea behind the sub- 
routines for mass and stiffness is outlined here for continuity. 

Mass and stiffness matrices of the complete structure , (fluid and 
bladder) are calculated using finite-element approach. In this approach, 
a continuous structure (fluid and bladder each separately) is assumed to 
be composed of simple, small structural elements - the so called finite 
elements - such as tetrahedron, pentahedron, triangular plates, quadri- 
lateral plates, etc. The procedure to obtain the finite-element mass and 
stiffness matrices is based on kinetic and strain energy principles, 
respectively. 
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The kinetic energy for a complete structure may be expressed as 

T= % fff • S (*» dx <*y dt ( 66 ) 

where: T = kinetic energy 

= mass density 

S> *= time rate change of deflection 

t “ time 

x,y,z - global coordinates 

The difficulty in integrating equation (66) is expressing the deflection 
£ (x,y,z,t) as a continuous function over the complete structure. In the 
finite-element approach, however, this apparent difficulty is circumvented 
by idealizing the structure to be comprised of many small structural ele- 
ments for which £ (x,y,z,t) can be expressed as a continuous function. 
Thus, the expression (66) is valid for each of the finite-element of the 
structure. Then the kinetic energy of the structure is the summation of 
the kinetic energies of each of the finite elements, that is, 


T =21 T i (67) 

where i refers to one particular finite element "i". 

The common junction of finite elements is denoted as panel points, nodes 
or joints. Now, however, the deflection 5 (x,y,z,t) is easily expressed 
as a simple function of the joint deflections. These element joint de- 
flections are then generalized coordinates or degrees of freedom of the 
complete structure. 

The approach is to derive the mass matrix for finite-element, n i", 
in a convenient local coordinate system and then transform it to the 
Global coordinate system. The technique is outlined here: 

iz=i{V^} T W x (68) 

A* 

where W; = the mass matrix in the local coordinate system for the 

ith element. This mass matrix is obtained by integra- 
tion using an assumed displacement function. The 
discussion is deferred till later. 
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- the time rate of change of the joint deflections of 
finite- element, "i". This is in local system. 

The deflections in the local coordinate systems are related to deflections 
in the global coordinate directions by a transformation matrix, fT*! * 

Thus 



= [Tim 


(69) 


where 


= the joint deflections of finite element, "i", i 
^ the global coordinate system. 


Using equation (69) in equation (68) 




where 


■ cn 


(70) 

(71) 


is the mass matrix with respect to the global coordinate system for the 
ith finite-element. Further, all the elemental mass matrices are finally 
assembled to give the mass matrix of the total structure, as shown in 
equation (67). 

The development of the finite- element stiffness matrices is similar 
to that of the mass matrices. The strain energy for the structure may 
be expressed as the summation of the strain energies of each finite 
elements. That is, 


U=2LUi 


(72) 


As was done for the finite- element mass matrix, the stiffness matrix for 
finite- element, "i", is derived in a convenient local coordinate system. 
Thus , 




( 73 ) 
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where J * the stiffness matrix with respect to local coordinate 

* directions for finite element, "i”. This stiffness matrix 
is obtained by integration using an assumed displacement 
function. This will be discussed later. 



the joint deflections of finite element, i, measured in 
local coordinate system. 


The same transformation matrix, cn , which was used in equation (70) 
is used here to relate the deflections in local coordinates to deflec- 
tions in global coordinates. Substitute then, 


% [h® fallal 

where -j" 


(74) 


(75) 


is the stiffness matrix with respect to the global coordinate system for 
the ith element. 

Euler angle rotations at some joints (where the body coordinate is 
needed to be different than that of the global coordinates) are input in 
the program to allow the joint degree of freedom at these points to be 
different than that of global x,y,z directions. 


However, in case of fluid, there is another item to be taken care 
of as far as stiffness matrix is concerned. This has to do with the 
surface elements of the fluid. The item concerned is known as gravita- 

O 

tional potential. The energy contributed is known as the gravitational 
potential. This is caused by fluid movement in the gravitational field. 
Development of the gravitational potential stiffness effect is given in 
Section 2.2.4. 
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2.2.3 Stiffness and Mass Matrices In Local Coordinates (Solid 
Element) 

Mass Matrix (Triangular Element) 



Triangle Plate Membrane Triangle Bending Triangle 

Element Plate Element Plate Element 

Figure - 2 

The consistent mass matrix for a triangle plate element are calculated 
as a combination of the two following elements. 

(i) consistent mass matrix for a membrane triangle plate element, 

(ii) consistent mass matrix for a bending triangle plate element. 

The elements of Ae mass matrix for the membrane triangle plate element 
represent the distributed mass properties of the triangle. These matrix 
elements are calculated by assuming a quadratic displacement field^) , 

The elements of the mass matrix for the bending triangle plate element 
represent the distributed mass properties of the triangle. These matrix 
elements are calculated by assuming a cubic displacement field*^. 

Stiffness Matrix: (triangular element) 

The stiffness matrix for the triangle plate element is calculated in 
the same manner as the mass matrix and the technique and displacement 
fields are exactly the same. 

Mass Matrix (quadrilateral element) 



Figure 3 
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This mass matrix is calculated by taking the average of the four over- 
lapping triangles created by the diagonals (1.3) and (2.4) « The triangles 
are handled as discussed before. Thus, the quadrilateral case is nothing 
but a combination of triangular case. 

Stiffness Matrix (quadrilateral element) 

This is handled in exactly same fashion as the mass matrix case. 

2.2.4 Stiffness and Hass Matrices in Local Coordinates (fluid 
finite element) 



for Tetrahedron Element 


The basic fluid element is a tetrahedron; pentahedron elements and 
hexahedron elements are synthesized, simply by placing six and ten 
overlapping tetrahedrons together respectively and averaging the result. 
The averaging is carried out to eliminate the bias, if any. 

For each tetrahedron element, a local cartesian coordinate system 
is defined so that vertex 1 is the origin, the x-axis includes vertex 2, 
vertex 3 lies in the x-y plane and vertex 4 always has a positive Z- 
coordinate (figure above). 

This element considers a linear displacement field (constant strain). 
This is boundary conformable. The displacement field throughout the 
element is expressed in terms of coordinate locations and appears as 

l 

Zr (x>y> £.t) — 5- 0 -+ 5# 


(76) 
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The coefficients *a k (t) , K = 0, 1, 2, 3 are eliminated in terms of the 12 
vertex displacements. 

The mass matrix for the fluid elements is obtained by expressing the 
kinetic energy 

T - ^2 f <=*v (77) 

Wat 

where: T = kinetic energy 

if = i 0 +3£c-f 

P = mass density 

This gives rise to a (12x12) mass matrix. 

The stiffness matrix for the fluid element is obtained by expressing 
volumetric dilatation strain energy and gravitational potential energy in 
terms of vertex displacement coordinates 



where 


Tjr - l/M f C 67 '^) C 67 * *) ** 09 ) 


u. 


D 


n 

"e 


= volumetric dilatation energy 
= gravitational potential energy 
= fluid bulk modulus 
- volumetric strain 
= unit outer normal 

= a unit vector parallel with the gravity vector g, 
opposite sense, i.e., e - -g/g 


but of 


An observation is made with respect to the gravitational potential energy 
as expressed in (79). Since it is a boundary conformable element, the 
surface integrals such as (79) will all cancel each other throughout the 
interior of the fluid in a container, since n on common element boundaries 
is equal and opposite. Thus, the gravitational potential energy will de- 
pend only on displacement coordinates at the boundary of the entire volume 
of the fluid, the free surface, the wetted container wall and also, in 
this case, the bladder. Also, for a rigid tank, VCf • "n is non- zero only at 
the free surface where *e * ri; thus. 
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- 14 (80) 
f rrc At(rr|'Ci 

Stiffness coefficients corresponding to gravitational potential 
and volumetric strain energies are thus derived. 

2.3 Mechanical Equivalent - The mechanical equivalent in this case has 
been calculated on the following basis. The forces and moments developed 
by the model due to an external disturbance should correspond to the 
forces and moments exerted by the fluid under similar conditions. 

Since any arbitrary liquid motion can be thought of as superposition 
of different slosh modes it suggests itself that the mechanical equi- . 
valent must consist of a series of spring mass systems, the masses 
corresponding to the effective amounts of liquid oscillating in different 
slosh modes. The frequencies of the mechanical equivalent must correspond 
to the frequencies of the elastomer for the mode they represent. However, 
to account for the part of the fluid that does not take part in the motion, 
one may also add a mass without a spring. 

A short definition of the mechanical equivalent is as follows: 

The equations of motion for a base driven elastic system are 

{ i) 

where : 

^ : modal coordinates, 

% : modal velocities, 

^ : modal accelerations, 

4 ’ ' discrete mode shapes, 

T : rigid body transformation matrix, 

% : discrete base accelera tions (6 degrees of freedom), 

^ : modal damping, 
c<> : frequencies, 

: mass matrix 
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The base reactions to support a given modal accelerations are 

J t. 7 

Noting that 


^ c< {+] T [n] [t] {it} 


( 82 ) 


(83) 


Substituting from(S3) in(82) 

T 


T 


{n hhki 

where the proportionality parameter is a function of frequency. 

The product 

"J" 

is a (6x6) matrix, in 

general, for a given mode. This may be called a matrix of forces and 
moments. The equivalent slosh parameters can be interpreted as shown 
in the following sketch. 


(84) 



where 

n^: effective slosh mass 

x: distance of line of action from the reference point. 

The reference point is the same as the reference for the rigid body 
transformation. 

This equivalent model is only valid for loads at the reference 
point. This may not be used for getting any detailed information 
concerning what is happening inside the elastic system. 



26 


In this particular case, the final product is a matrix (3x3) . The 
term (2,2) is the slosh mass in the y-direction and the ratio of terms 
(2,3) and (2,2) is the distance x. It should be noted that x may be 
outside the tank boundary. The slosh mass and x are printed out in the 
computer run. The slosh stiffness is the frequency squared because the 
modal displacements are normalized such that the generalized mass is 
unity. To obtain the mechanical equivalent slosh parameters of the pre- 
ceding sketch, only y- translation and 0 Z rotation need be considered in 
the rigid body transformation matrix T . Thus, the general 6x6 slosh 
matrix is a 2x2 matrix for the particular case considered here. The 
1,1 term is the slosh mass in the lateral y-direction. and the ratio of 
terms (1,2) and (1.1) is the distance x. It should be noted that x may 


be outside the tank boundary. The slosh stiffness is the slosh frequency 
squared time the slosh mass, m . 
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3. MODEL 

The modeling -of the system is achieved in two distinct parts. First, 
the static free surface of the fluid/bladder system, within the constraints, 
is established (refer to 2.1). This static free surface is subsequently 
used in modeling the total system as an assemblage of finite-elements for 
vibration. 

3 *1 Finite .Element Model of the System - This is an axi symmetric sys- 
tem. The system thus can be subdivided into four quadrants., 1, 2, 3 and 4, 
as shown in the figure. Quadrants 2, 3 and 4 are reflected on to quadrant 
1 with the following technique whereby the analysis of only one quadrant is 
necessary for the dynamic analysis of the total system. This renders the 
problem relatively small as far as sizes of the matrices are concerned. 



Figure 5. Quadrants and Coordinates Systems 
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x y z is the global coordinate system (right handed) 

x v Zr is the so called body coordinate system (right handed) 

& J h b 

The global and body coordinate system systems are t o different coordinate 
systems related through the Euler rotations. The global system remains 
fixed but the body coordinate system varies for different points. 

The global system is used for interior fluid and the body coordinate sys- 
tem is used on the boundary. The use of body coordinate system on the 
boundaries makes the handling of the boundary condition easier. The body 
coordinates are shown in the figure for a few important points. The re- 
flections of quadrants and the derivation of the Euler angles are discussed 
in the following pages. 

• 3.2 Reflections of Quadrants 2 . 3 and 4 on Quadrant 1 for 
Lateral Slosh - For the case of lateral slosh and the excitation direction . 



Figure 6 
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■where 

and 


Thus , 


v and w are the velocities in the global x, y and z directions 
a typical point in quadrant 1 (on the boundary) 

a typical point in quadrant 2 (on the boundary) 

a typical point in quadrant 3 (on the boundary) 

a typical point in quadrant 4 (on the boundary) 


q’s: local coordinates 

A: a typical point in quadrant 1 (interior) 
B: a typical point in quadrant 2 (interior) 
C: a typical point in quadrant 3 (interior) 
D: a typical point in quadrant 4 (interior) 
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( 86 ) 


(87) 


it is apparent 


- K 


( 88 ) 

( 89 ) 
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There exists a relation between <5* s and global coordinates. They are 
(writing them in terms of their velocity or displacement components) 


fll 

f U 'j 


J v I 


i C*) 

l* 

J ** 1 

\s 

J 

U»J 


(90) 




(91) 


Recognizing the relations (88) and (89), the other quadrants relations 
are not followed through any further. From relations (85) and (91) the 
relation between the degrees of freedom of quadrants 1 and 2 for interior 
of the fluid can be written as 


Thus 



zl 

' 









n 






h— 

□ 














"1 






v 

u 

Or 


W*- On] {*} 


( 92 ) 


( 93 ) 


where 



[ T « 1 
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If the boundary degrees of freedom were also in the global comparison 
system, the same transformation could be applied to them. But as the 
situation stands, they are in the body coordinates system which in 
general will have a different transformation than [ T 2 il • 
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4. RESULTS 

The results of the study has been obtained in three distinct steps. 
They are: (1) Static Equilibrium Shape, (2) Vibration Analysis and, 

(3) Equivalent Slosh Mass for Dynamic Characteristics of the System. 

The static equilibrium shapes obtained in the first part are fed into 
the second part for establishing the grid pattern for finite-element 
analysis of the system. 

4.1 Static Equilibrium Shapes - The technique used for this has 
been discussed in Section 2.1. A computer program has been developed 
for this purpose which is listed in Section 7. 

The static equilibrium shapes for fill conditions of 80%, 60%, 50%, 
40% and 20%> for lg case have been obtained using the available computer 
program. However, for the case of 10"^g conditions, only fill ratios of 
80% and 60% have been successfully run through the program. The other 
cases of 50%, 40% and 20% were extrapolated from available results. 

This was achieved with the help of graphical representation between the 
percentage fill and the coefficients of the power series, assumed for 
the static surface, see Figures 11G through 14G. The graphical repre- 
sentations of the static free surface and the coefficients are presented 
in Figures 1G through 14G. 

The values of other parameters used in all of these are: 

(i) Radius of the sphere = 8 inches (20.32 cm) 

(ii) Young's modulus for the bladder material = 200 lb/in 2 

(1.37895lxl0 6 N/m 2 ) 

(iii) Thickness of the bladder = 0.06 inches (0.1524 cm) 

(iv) Mass density for mercury = 0.0013 lb/in/sec 2 (0.00059 Kgm) 

(v) Ullage pressure =0.10 lb/in^ (689.48 N/m 2 ) 
and the assumptions made are: 

(vi) Bladder attachment level - diametral plane 

(vii) Order of the polynomial = 4 

Figures 11G through 14G may be used to establish the coefficients 
which define the static free surface for any gravity and any percentage 
fill without making any more computer runs. 
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x 



STATIC FREE SURFACE FOR 80% FULL, lg 
FIGURE 1G 




STATIC FREE SURFACE FOR 60% FUU-, lg 


FIGURE 2G 
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x . 



STATIC FREE SURFACE FOR 40% FULL, lg 


FIGURE 4G 
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* 



STATIC FREE SURFACE FOR 20 % FULL, lg 


FIGURE 5G 



IGURE 6G 
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% 



STATIC FREE SURFACE FOR 60% FULL, 10 _5 g 


FIGURE 7G 



STATIC FREE SURFACE FOR 50% FULL, 10 _5 g 


FIGURE 8G 
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STATIC FREE SURFACE FOR 20 % FULL, 10~ 5 


FIGURE 10G 
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4.2 Vibration Analysis - The iterative Rayleigh-Ritz technique used 
to calculate the vibration mode shapes and frequencies has been described 
in Section 2.2.1. The computer subroutine m TM0DE2 u used to calculate these 
modal values already exists at the MSFC computer center and is, thus, omitted 
from the program listings of Appendix A2. 

The calculated vibration mode shapes were used in the calculation of 
the mechanical equivalent slosh parameters which are given in Section 4.3. 

The calculated vibration frequencies are listed in Section 4.3. 

4.3 Mechanical Equivalent Slosh - The technique used to calculate the 
mechanical equivalent slosh parameters has been described in Section 2.3. 

The computer subroutine tf MECHEQ" used to calculate these parameters is in- 
cluded in Appendix A2 . 

A great deal of difficulty was encountered in the identification of 
the slosh modes among the many vibration modes calculated. It was decided 
to identify as slosh modes those modes having the largest slosh mass. The 
following tables give the slosh frequencies, masses, stiffnesses and attach 
stations. It should be noted that numerical precision of the computer 
casts some doubt on the accuracy of these modal results. 
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Table 4.3-1 

SLOSH PARAMETERS 
80% FILL CONDITION 
lg ACCELERATION 


f 

(Hz) 

Mass 

lb-sec^/in 

Stiffness 
lb /in 

Attach Station (ic) 
in 

1.0756 

.0037 

.1690 

4.2 

1.0828 

.0022 

.1018 

10.0 

1.0946 

.0002 

.0104 

6.8 

1.1193 

.0001 

.0044 

8.9 

1.1235 

.0014 

.0698 

5.0 

1.1337 

.0011 

.0558 

6.2 

1.1709 

.0008 

.0411 

1.5 

1.2429 

.0180 

1.0978 

6.0 

1.2471 

.0014 

.0860 

4.3 

1.2640 

.0250 

1.5768 

3.9 

1.3188 

.0130 

.8926 

5.7 

1.3276 

.0035 

.2435 

4.8 

1.3709 

.0120 

.8903 

5.1 

1.4495 

.0007 

.0581 

11.4 

1.4625 

.0003 

.0262 

8.4 

1.4729 

.0350 

2.9976 

3.7 

1.4952 

.0001 

.0097 

-4.0 


f 

Table 4.3-2 

SLOSH PARAMETERS 
60% FILL CONDITION 
lg ACCELERATION 

Mass Stiffness 

Attach Station (x) 

(Hz) 

lb- sec in 

lb/ in 

in 

1.0805 

.0008 

.0350 

22.7 

1.1020 

.0093 

.4458 

. 2.0 

1.1284 

.0660 

3.3176 

5.1 

1.1384 

.0002 

.0112 

7.3 

1.1741 

.0230 

1.2516 

7.7 

1.1871 

.0007 

.0384 

-3.8 

1.2307 

.0041 

.2452 

3.7 

1.2355 

'.0550 

3.3143 

3.5 

1.2412 

.0240 

1.4596 

1.7 

1.2647 

.0830 

5.2406 

1.5 

1.2748 

.0056 

.3592 

-1.3 

1.2764 

.0420 

2.7014 

0.8 

1.3093 

.1400 

9.4748 

2.3 

1.4149 

.0000 

.0001 

-94.0 
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Table 4.3-3 

SLOSH PARAMETERS 
50% FILL CONDITION 
lg ACCELERATION 


f 

Mass 

Stiffness 

Attach Station (x) 

ISsl 

lb-sec2/ in 

lb/ in 

in 

1.5270 

.0003 

.0285 

-3.6 

1.5692 

.0062 

.6027 

10.2 

1.5888 

.0060 

.5979 

12.3 

1.8572 

.0002 

.0245 

-2.6 

1.9030 

.0180 

2.5735 

9.6 

2.0489 

.0160 

2.6518 

6.7 

2.0725 

.0000 

.0100 

4.7 


f 

Table 4.3-4 

SLOSH PARAMETERS 
40% FILL CONDITION 
lg ACCELERATION 

Mass Stiffness 

Attach Station (x) 

(Hz) 

lb- sec2/ in 

lb/in 

in 

1.0594 

.0004 

.0199 

3.3 

1.0988 

.0001 

.0028 

-3.3 

1.1035 

.0210 

1.0258 

1.3 

1.1071 

.0008 

.0397 

5.0 

1.1674 

.0004 

.0200 

-2.8 

1.2077 

.0017 

.0978 

10.4 

1.2184 

.0001 

.0068 

35.6 

1.2429 

.0002 

.0127 

-0.3 

1.2570 

.0005 

.0336 

28.0 

1.2705 

.0001 

.0072 

25.0 

A jump was 

made to higher frequencies 


5.4081 

.0003 

.3938 

37.2 

6.0773 

.0000 

.0032 

11.4 

7.0748 

.0000 

.0000 

-177.0 

8.1944 

.0004 

1.0516 

28.3 

8.2019 

.0001 

.2094 

27.5 
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Table 4.3-5 


SLOSH PARAMETERS 
20% FILL CONDITION 
lg ACCELERATION 


f 

Mass 

Stiffness 

Attach Station 

iSsl 

lb-sec2/ In 

lb/ in 

in 

1.1170 

.0019 

.0936 

3.7 

1.1209 

.0060 

.2971 

3.1 

1.1250 

.0017 

.0084 

-0.4 

2.6893 

.0099 

2.8272 

1.9 

3.3176 

.0000 

.0162 

-16.6 
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Table 4.3-6 . 

SLOSH PARAMETERS 
80% FILL CONDITION 
10- 5g ACCELERATION 


f 

Mass 

Stiffness 

Attach Station Cx) 

(Hz) 

lb-sec^/in lb/in 

in 

0.9374 

.0390 

1.3469 

3.7 

0.9803 

.0055 

.2080 

1.0 

1.1229 

.0006 

.0319 

17.0 

1.1420 

.0022 

.1151 

-2.8 

1.1697 

.0023 

.1235 

9.6 

1.2470 

.0000 

.0011 

22.0 

1.2914 

.0012 

.0823 

9.4 

1.3485 

.6025 

43.2559 

7.7 

1.3604 

.3385 

24.7301 

7.3 

1.3964 

.0113 

.8699 

9.2 

1.3967 

.0064 

.4940 

3.4 

1.4706 

.0011 

.0981 

2.2 

1.5501 

.0037 

.3546 

4.6 

1.6660 

.0003 

.0285 

Table 4.3-7 

SLOSH PARAMETERS 
60% FILL CONDITION 
10- 5g ACCELERATION 

-26.2 


f 

Mass 

Stiffness 

Attach Station Cx) 

(Hz) 

lb-sec^/in 

lb/ in 

in 

0.9789 

.0006 

.0246 

18.2 

1.0190 

.0033 

.1334 

. 5.4 

1.0461 

.0530 

2.2743 

6.4 

1.1590 

.0040 

.2098 

5.9 

1.2204 

.0001 

.0072 

8.1 
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Table 4.3-8 

SLOSH PARAMETERS 
50% FILL CONDITION 
10“ 5 g ACCELERATION 


f 

Mass 

Stiffness 

Attach Station (H) 

(Hz) 

lb-sec2/in Ib/in 

in 

1.0629 

.0008 

.0380 

-10.9 

1.1177 

.0180 

.8636 

-0.4 

1.1392 

.0043 

.2197 

-9.9 

1.1759 

.0095 

.5195 

0.7 

1.1766 

.0030 

.1620 

-2.8 

1.2032 

.0094 

.5393 

0.3 

1.2247 

.0000 

.0002 

17.9 

1.2859 

.0030 

.1972 

-2.4 

1.3917 

.0000 

.0039 

5.0 

1.3942 

.0020 

.1569 

-2.6 

1.3949 

.0027 

.2103 

3.8 

1.4356 ' 

.0030 

.2431 

-2.0 

1.5645 

.0003 

.0320 

-10.6 

1.5788 

.0022 

,2118 

2.8 

1.6007 

.0001 

.0066 

-22.0 



Table 4.3-9 




SLOSH PARAMETERS 




40% FILL CONDITION 
10“ -5 g ACCELERATION 


f 

(Hz) 

Mass Stiffness 

lb-sec^/in Ib/in 

Attach Station Cx) 
in 

1.0036 

.0038 

.1519 

-6.3 

1.0191 

.0001 

.0021 

-39.0 

1.0238 

.0006 

.0212 

2.9 

1.0384 

.0026 

.1099 

-8.7 

1.0660 

.0015 

.0685 

1.4 

1.0973 

.0010 

.0475 

17.7 

1.1022 

.0036 

.1745 

13.3 

1.1441 

.0001 

.0078 

42.9 

1.1496 

.0004 

.0200 

12.1 

1.1772 

.0045 

.2457 

2.9 

1.1823 

.0220 

1.2059 

9.4 

1.2267 

.0000 

.0001 

10.3 

1.2270 

.0008 

.0471 

0.7 

1.2514 

.0250 

1.5727 

3.6 

1.2555 

.0027 

.1698 

20.1 

1.2879 

.0060 

.3929 

4.7 

1.3154 

.0008 

.0548 

-0.8 

1.3403 

.0001 

.0070 

-3.5 

1.3754 

.0004 

.0306 

-0.5 
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Table 4.3-10 

SLOSH PARAMETERS 
20% FILL CONDITION 
10-5g ACCELERATION 


f 

Mass 

Stiffness 

Attach Station Csl 

ISal 

lb- sec/in 

lb/ in 

in 

1.0284 

.0008 

.0314 

-4.9 

1.0661 

.0030 

.1363 

-2.6 

1.1009 

.0000 

.0003 

-2.07 

1.1260 

.0000 

.0006 

55.9 

1.1440 

.0012 

.0615 

2.0 

1.1624 

.0003 

.0154 

23.5 

1.1703 

.0017 

.0939 

1.1 

1.1933 

.0004 

.0235 

-9.7 

1.2512 

.0000 

.0025 

-48.5 

1.3078 

.0000 

.0023 

7.5 

1.3175 

.0052 

.3592 

-1.5 

1.3381 

.0019 

.1320 

5.2 

1.3401 

.0011 

.0746 

-10.1 

1.3642 

.0007 

.0510 

-29.8 

1.4709 

.0001 

.0048 

-67.8 
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5. CONCLUSIONS AND COMMENTS' 

It is felt that a significant contribution was made in this study 
to the state of the art in finite element fluid analysis. As with all 
new investigative analytical studies, review of the work performed re- 
veals that although significant advances were made, several items that 
should have been studied further or have been performed in a different 
manner if time had permitted. 

Several different approaches were developed and programmed to cal- 
culate the static equilibrium shape of the fluid/bladder system. Most 
of them did not perform satisfactorily due to the numerical or convergence 
problems. Finally, being limited by the time available, a two-dimensional 
representation - an infinite channel - was chosen for the representation 
of the equilibrium shape as the only approach showing a reasonably good 
convergence. A future effort should review this approach critically and 
perhaps extend it to a 3-dimensional representation. 

A second item that should be investigated is the use of double pre- 
cision in the computer programs, particularly for the calculation of 
vibration modal properties. This study pointed up the possible need for 
double precision because a bladder with small modulus of elasticity was 
used in the analysis. Thus, clear separation of the fluid slosh modes 
from fluid circulation modes was clouded. In addition, use of low shift 
values, , (see Section 2.2.1) sometimes resulted in failure to de- 
compose the dynamical matrix ^ due to singularity. This is 

obviously a computer accuracy problem. 

The data generator computer subroutine, used to calculate joint X, 

Y, Z locations, degree of freedom values, Euler angles, and finite element 
joint numbers has some limitations which should be removed in future 
studies. One of these limitations is the requirement for vertical radial 
cuts. This was used to minimize user input but it became obvious later 
that odd shaped elements resulted. A more general data generator should 
be coded to allow more user control on the shape of the elements. 
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The representation of the lateral slosh by a spring mass system has 
been achieved on the assumption that the modes are completely uncoupled. 
A detailed investigation of slosh equivalent modeling techniques should 
be persued in this direction in a follow-on effort. 
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7 . APPENDIX 

The analysis of the system is accomplished in two main steps. The 
first is the static free surface definition and the second is the dynamic 
analysis. The dynamic analysis program includes the automatic generation 
of joint X, Y, Z values, degree of freedom numbers, Euler angles, element 
joint numbers, calculation of mode shapes, and frequencies, mechanical 
equivalent slosh mass and plots. A computer program has been developed 
for these steps. A schematic flow chart is shown here for the system 
analysis steps. The listing of the static free surface computer program 
is given in Appendix A- 1 and that of Dynamic analysis computer program in 
Appendix A-2. The important parameters input to the programs and subrou- 
tines are explained in Appendices B-l and B-2, respectively, using typical 
input listings. 

A brief summary of important subroutine functions are presented in 


the following pages. 



Brief Summar y: 



o> 

o 
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Summary of Programs : 

STATIC FREE SURFACE obtains the static free surface 

DYNAMIC ANALYSIS obtains the dynamic characteristics of the 

system (frequencies of mode shapes) : and the 

mechanical equivalent 

Summary of Subroutines (used in the DYNAMIC ANALYSIS Program) : 

GNXYZ3 ■' input data for 'FINEL 1 

XYZEUL automatic generation of input for 'FINEL 1 

* 

FINEL generates mass and stiffness matrices 

FLUID generates mass and stiffness for fluid only 

GRAVTY generates gravity contribution to stiffness matrix 

TRNGL generates mass and stiffness for non-fluid triangular 

elements 

QUAD generates mass and stiffness for non- fluid quadri- 

lateral elements 

MODES obtains frequencies and mode shapes 

MECHEQ obtains mechanical equivalent for the sloshing fluid 

and bladder 


PLOT 


plots the mode shapes for the mid- plane 
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C- 

C- 
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2001 
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c- 
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c- 

c- 
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c- 

c- 

c 
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c 
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APPENDIX - A1 

STATIC FREE SURFACE PROGRAM 


PROGRAM StPS ( INPUTtOUTPUT*TAPe'S=INPUr* rAPt6*0UrPUn 004840 

IMPLICIT DOUBLE PRECISION CA-H*0-Z) 004850 

004860 

DIMENSION F(15), DFOQ ( 15, lb) • C<i6,16)* khSVUd), VF.C(lb)* 004870 

1 VEC'i(lS), €11(15,15), CI2U5,15)§ wORK<ib*15), 004860 

2 T ( i 5 , 1 5 ) , rT(15*15), ClC2(lbtlb)',PLriND(83) *PLTDEP(83) 004890 

004900 

EQUIVALENCE (PLIINO ( 1) ,WORK U ) ) , (PLTULP(l) , WUHK (84) ) 004910 

004920 

EVTENNAL f T1 ,FT2»F I 3 ,ET4,FT5>E T6,f T 7 004930 

0049^0 

COMMQN/BLK1 / E,Th ,KH0,6,R*K *jmO 0 04950 

COMMON U ( 1 1> > 004960 

004970 

FORMAT f 15,301 /> B) 004980 

FORMAT ( lrtl , 15 (/) ,2 dX,*maxTMUH l TER AT I ON LIMIT REACHED*) 

FORMAT ( IX, 5D2 0.10) 

005000 

DATA K Q / 15/ 005010 

nATA A/0. 0000/, 8/6.0000/, FPS/l.OD-J/, E/20U.0DQ0/ 

DATA TH/O .060 00 / , RMO/0 . OO i Ju 00/ » 3/8.01)00/ 

O A T A I N/ 1 ft / , N‘C/4/ 005040 

005050 

V AR I ARLES- 

IN = NO. OF INTERVALS F 0 K NUHFhICAL I N TEGR 4 T 1 ON . 005060 

A = LOWER LIMIT 005070 

9 ■= UPPER LI'UT 005080 

MO = ORDER Of THE POLYNOMIAL 0 05090 

TH r THICKNESS OK I HE R L A 0 L) t K • 005100 

F - YOUNGS MODULUS FOR THE H L A 0 D L R MATERIAL, 005110 

rho = MASS DENSITY OF THE FLUID. 005120 

G = ACCELERATION out GRAVITY 005130 

P = ULLAGE PRESSURE. - 005140 

R = RADIUS OF THE RARRtr't. 005150 

NT = WO. OF CONS TRAN l EQUATIONS x 005160 

IJV = ULLAGE VOLUME 005170 

FP5 S EPSILON TO CUMpAKF OU(I>-$. 005180 


$■ fr ■» •» ■» * v- ■» ■> a * * «• & a a v ■» # ->4- # # «■ ■» a * -a * # » ■» * •» -a- a k -jj * «• * -it * * a ■& -a- * * # * » * * * #. * * * * «■ a 


INPUTS- 

MO, G,P, PC T (15,3017.8) 


005190 

CAM START 005200 

P F A O ( 5 » 1 V 0 ] ) NO , G * P , PC T 005210 

Call 7Fku ( DFi'.'Q , KO »KN»KO ) 005220 

CAM. 7ER > (C , KO , 7\<ji<Q ) 005230 

C Ai L /FRO ( Vr. C , a w , 1 , K >) ) 0 0524 0 

p r = 2 . 0 * A 3 I N ( i . o ) 005250 



n n n non n o o no 
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APCVAL = 1,00 * PI * R / 2.U 

IJV = 0,5* ( 1 • O-Pcn 
vfc (mo = uv 

Nl ss NO + 1 
N? = NO ♦ ? 

GFNFRaTE The MATRIX (C)0. 

FTPST HO* 

C ( 1 * 1 ) = 1 • o 

no 130 J=1*N0 
130 CdiJ+1) = B**.J 
SFCOMD HO* 

C<?«1 )= o.u 

c < ? • 2 > * l . a 

no 135 J=?*NO 

135 C(?*J+1) = B**(J-l)*FLOflr(J) 

TNI PD HO* 

C ( 3 , 1 ) = 0 • 0 
C < 3*21=1.0 
no mo j=3*»u 
i4n cn»j) =o.o 
fourth ROW. 
no 145 J = 1 * N X 

145 C ( 4 . J ) = / FLOAT (J) 

CALL WRITE (C*NC*N1 *2HCO«<U) 

PART IT TON (C) TO MAKE (CU) 4HD (Cl 2) 

nn iso I = UMC 
nn ]50 J=1 *nC 
l50Cl!(IfJ)=C(I*J) 

MCI = NiC + 1 
NJ) MC = N1-NC 
no 155 1=1 *NC 
on 155 J = 1 * N 1 N c 
155 Cl ? ( I* J) = C< I * J + NC) 

CALL WrtPF (C i 1 *nC*NC* 3HC1 1 *KO) 

CALL WRITE <Cl2*NC*NlNCf3HC12,KQ) 

GFT (Cll) INVERSE 
CALL TNV 1 (C n • WORK f iMC«KO ) 

CALL WHITE ( WORK* NC*NC»6HC1 1 lNVtR-0) 

GFT {-Cl 1 ) INVERSE times ( C 1 2 ) 

CAi L MULT { work fC12* ClC2f NC *NC *NiNC ♦ KG #KQ) 

CALL ZERO ( T «KQ *KQ »KQ1 

nn 160 1 = 1 » 4 

nn 160 J= 1 « N 1 NC 
160 T ( T , J) =-ClC2 ( 1 , J) 
no 165 1 = N C 1 » N 1 
TMC!=I-NC 
165 T (T *INCi ) =1 . 0 

CALI wRITF ( T*Nl *NlNCt4HTR'ANfK(j) 
c~ FORM MATRIX (Cll) INVF.PSF Tints THE RIGHT HAND SI OF 


005290 

005300 

005310 

005320 

005330 

005340 

005350 

005360 

005370 

005380 

005390 

005400 

005410 

005420 

005430 

005440 

005450 

005460 

005470 

005480 

005500 

005510 

005520 

005530 

005540 

005550 

005560 

005570 

005580 

005590 

005600 

005610 

005620 

005630 

005640 

005650 

005660 

005670 

005680 

005690 

005700 

005710 

005720 

005730 

005740 

005750 

005760 

005770 

005780 
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C&Ll 7ER0 { Q » Kij » i *K9) 

005790 


CM. L v<ULV ( WORK , V£C • NC* NC » 1 »KU»Ku) 

005800 


0 t N? > - S I MP S ( A ♦ ft * IN , F T 7 * T UUM , I HUM ) / ABC V At - 1 . 0 

005810 


C A 1 L WRI TE (U*N2* 1 *6HO ( I ) -I iK{jj 

005B20 

c 


005830 

c 


005840 

c- 

LOOP TO DETERS TNti (J(I)-S. 

005850 


MM - 0 

005860 

100 

CONTINUE 

005870 


Cfll.L MULT (C*«*«0HK*NC*N1 *1 *K‘J*KO)- 

005880 


C*( L WRITE ( WORK « NC f 1 f 6HC0 NS TN,kQ) 

005890 


WRITE <6* 2000 ) MM 

005900 

2000 

FORMAT (* ITERATION NO. *,73) 

005910 


MM s MM +1 

005920 


T f (MM.GT.50) 60 fO 190 


c- 

FOP thf first n*i«f-s. 

005940 


no no L = 1 ♦ n l 

005950 

110 

F ( ! ) = SIMPS < A ,B, IN,F T1 ,L* IOUM) 

005960 

c- 

FOP THE LAST N + 2 N D « F 

005970 


C A | L WRITE It *Ni * 1 ,ZHFN,K<)> 

005960 

c- 

CM CULATE (DFOj) matrix. 

005990 


do ] 1 5 1 = 1 , N 1 

006000 


DO 1 IS J=1 , 1 

006010 

1 1 s 

PFOOd.J) s SlMf>3(A*fi,INI*fT3M,J) 

006020 


on 120 I = 1 » Ns 1 

006030 


on )20 J S 1 » l 

006040 

120 

OFOO(Jd) = OrOO(hJ) 

006050 


CHL WRITE ( OF DO » N 1 , N 1 « 4H0 E ON , * q ) 

006060 

c 

> 

006070 

c- 

OFT (T) TRANSPOSE AND MULTIPLY 

006060 


DO 1 7ft 1=1, Ml 

006090 


00 170 J=1,N1NC 

006100 

170 

T T ( J , T ) = f(hJ) 

006110 

c 


006120 


Cat L MULTR {TT,EtNlNC«Nl,l,Kvi,KQJ 

006130. 

c~ 

■OBTAIN' < T ) TRANSPOSE X (OFDU) X <T) 

006140 


CAU B T Art A <L>FOQ*T,Nl«NlNC'KU,K0) 

006150 


CM L WRITE (OEf)(jtNlNC*NlNC*EiHTTfJF(jT,Ky) 

006160 


C At | I NV 1 ( OFOU , T T , N 1 NC , K(J) 

006170 


cat L MIJL f 3 (IT , F » hi 1 N C * N I M C * t , * u . K Q ) 

006160 

c- 

OBTAIN 0EITA-A2 

006190 


C A | | , WRITE ( E * N 1 N C * 1 ♦ rt m 0 0 2 * K w > 

006200 

c 


006210 

c- 

C?FT delta- a. 

006220 

c 


006230 


CAM. MULT ( 7 ♦ F , }j w , N 1 , Ml N C , 1 * r; < .j , K w ) 

006240 

c 


006250 


C A | | w R I f E ( I/O i N 1 * i » 2 Hr' 0 , <* N ) 

0062b0 

c 


006270 

o 

GFT Q ( PRESENT ) = 0 (PREVIOUS) +00 

006280 

c 

. 

006290 


00 1 7S I - i • N 1 

006300 



non ooonooon nooo 
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17=, Q ( T ) = QU> - HQ ( I ) 

Q(N?) = SI^PS(A*tt* lN f FT7* IL)U»’1* JDUrt) /AHCVAL - 1,0 
C 

CAJ t WRITE OJ»N 2 « 1 *ShNOI-S*KU) 

W R T T fc" (6*2002) <U ( I ) $ 1 = 1 »N2) 

C 

on 180 Ltsl*Nl 

IF (OA 0 S(OQ<LL) ) •GT.F.PS) 00 To 100 
1 P 0 CONTINUE 

C- pi OT free SURFACE shape 
PLTIND (I) = R 

PI.TOEP(I) a 2 *P 
PI. T I NO ( 33 ) = 0 , 1 ; 

PLT0EP(B3) = 0.0 
DEL = R / 80.0 
X = A - UHL 
DO 188 1 = 2*81 

X r X ♦ OFL 
PI TINH ( I > = X 

PLTDFP(I) = Z ( X) 

188 CONTINUE 

PI.. T T NO { 82 ) = '» 

PI TDF.P(82) = ^ ( R ) 

CA|...L WRITF (PLTlNinb3,2,6HSHARc ,oj) 
no TO 195 
190 WPTTE (6*200 1 ) 

198 CONTINUE 
no ro 1 0 

FMO 

DOi JRL EPREC I S I ON FUNCTION SI HRS. ( 4 , H * N * F * NP * I NT ) 

IMPLICIT double PRECISION <A-h,u-2> 

THIS FUNCTION SUPPLIES 'ME DUMP W I CALL Y INTEGRATED VALUE OF THE 
jmtegrano. 

COmmON/RLK 1 / E * T H * R h 0 * G * P * P * NO 

COMMON 0(16) 


A -I. Owfr R LIMIT 
8—UPPE R LIMIT 
M-NO . OF INTERVALS 
f-funct ION 

HP-GENERAL I2FO COORDINATE NUMBER (as TO WHICH one IT IS) 

INITIALISE PARAMETERS 
TWnHa(9-A)/N 
M - TWOH/2.0 

sumeno=o. 0 

SU MM T 0 = 0 e 0 


TWOH- INTERVAL 
H-MALF INTERVAL 


006310 

006320 

006330 

006340 

006350 

006360 

006370 

006380 

006390 

006400 

006410 

006420 

006430 

006450 

006470 

006480 

006490 

006500 


006510 
006520 
Q06530 
0065*0 
006550 
006560 
006570 
006571 
006580 
006590 
006600 
00661 0 
006620 
006630 
006640 
006650 
006660 
006670 
006680 
006690 
006700 
006710 
006720 
006730 
006740 
006750 
0 Q67b0 
006770 
006780 



oooo o on ' noon 
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C 

C 


C. 


c 


c 

c 

c 


c 


SUwmEno-SUM of F(XI), I BEING 'EVEN EXCEPTING I-2N 006790 

SUWMID-S'JM OF F(*l>* I BEING ODD* 006800 

006810 

EVALUATE SUMEND AND SUMMIO. 006820 

On 1 K = 1 * N 006830 

X=A+FLOAT (K-l ) *rwOH , 006840 

SUvFNP = SUMENfi + F ( X *NP * INT) 006050 

1 SUmmjo = SUW.) ♦ F<X*H*NP* InT) 006860 

006070 

RETURN ESTIMATED VALUE OF THE INTEGRA T. 006080 

STMPS a ( 2. 0 *SIJmEnD + 4, n*SUMMln-F t A »NP» INT) + F (rt*NP» INT) ) *H/3. 0 006090 

006900 

RETURN 006910 

END 006920 

noiJRLFPRtClSlOU FUNCTION FT1 (X,nP,100) 006930 

IMPLICIT DOUBLE PRECISION {A-H,(W) 006931 

ThTS FUNCTION DEFINES THE INTEGRAND FOR F-S FUR DISCRETE VALUES OF006940 
FOP DIFFERENT EQUATION 'f-jUMBtRS CORRESPONDING TO DIFFERENT Q-S OR A00695G 
NP = EQUATION NUMBER 006960 

006970 

COWMON/8LK1/ EtTH«RHO«G*PvK*NO 006980 


COMMON vJ ( 1 6 ) 

NP =5 NO + ? 

T F (X ,Eu. 0.0 .UR. np . L E - 2 ) GO to b 

TFI..OT 1 = -FLOAT { (NP -1 ) »U\|P- 2 >) * ***<NP- 3 ) 
no TO 10 
6 TFi OT1 = 0.0 

10 CONTINUE 

IF (X ,Eu. 0.0 • uR « NP , LE • i) GO TO lb 

TFf. OT? = -FLOAT (NP-1) * x**(**P-2> 

GO TO 20 
IS TF( DTP = 0.0 
20 CONTINUE 

TFPM1= ( (t.*TH**3) / 12. 0) * ( (1.0 +ZP ( * ) *■*<?) **<-3. DO) ) 

1 * ( ( ( 1 . 0 + ZP < a) *«?) » (/OP { X) ) * TFLOU) 

2 ~H 2 .S 0 )*< 2 DP(X>* # ?) >* f F L U T 2 ) ) 

irpMp= (E*Ih*u ( n2 ) **2*0 # S>*( (( 1 # 0«-2P (X) **£)** (-0*5). ) 
1 * (ZP ( X) ) * TFLUT2) 

IF (X , F J , 0 , 0 ) GO TO 2 S 

jppM3= (HHO*G> ( X) * (- ( X** (N y -1 ) ) > 

jrpM 4 a P* ( X** (NP- 1 ) ) 

GO TO 30 
?S CONTINUE 

TFPM.3 = 0.0 
T PyM4 - O.o 


ORIGINAL PA G I,’ ro 

® «»* wmn 


006990 
007000 
007010 
007020 
007030 
007040 
007060 
007060 
007070 
007060 
007090 
007100 
0071 10 
007120 
007130 
007140 
007150 
007160 
007170 
007160 
007190 
007200 
007210 
007220 
007230 
007240 
007260 
007260 
007270 
007280 
007290 
007300 



no o non noon on 


AX- 6 


C 


C 


c 


3n FT1=TERM1MERM2MERM3 + TF:RM4 0 0 7310 

007320 

RFTUPN 007330 

FMD 007340 

DOtJRLFPKtClSION FUNCTION FT2 ( X , NP * I OU ) 007350 

IMO|. IC1T noUHlE PRECISION (A-H*(W) 007351 


THIS FUNCTION DEFINES THt INTECkaNO FOR DISCRETE X FOR (NO+2) EQUA007360 

007370 


COMMON/QLK 1 / E ■» T h * RhO ♦ 0 * R » R » NO 

COMMON 0 ( 1 S ) 

KIP = NO + ? 

TFPM* (F* fH*U (N<?) >*U1.0 + /P(X)**2>**(0.5) > 

FT? = TERM 
RETURN 

fMH 

DOUMLFPRECISION function Z(X) 

IMP| ICII OOURLF PRECISION (A-h*0-Z) 

THIS FUNCTION DEFINES THE FUNCTION Z(*> AT lU^CRE TE X. 

COMMON/3LK I / F 9 TH t RHO « G * ° i R * NO 
COMMON 0 ( 1 6 ) 

IF (R 4 t(j* X) GO TO 6 
7 s (R*#2-X**2) *5+R-W < X) 

CO TO 10 
6 7 - P - MX) 

in continue 

Rr TURN 
END 

OOUPLEPPECISICN FUNCTION MX) 

IMPLICIT double; PRECISION (A-h*(W) 

TH T S FUNCTION DEFINES MX) A I DISCRETE VALUE OF X. 

COMMON /BLK 1/ E * T H * RHO * G t P » R * «VO 
COMM O N O { 1 'o ) 

M] = NO ♦ 1 

W1 = o.o 

IF (X .EU. 0.0) GO TO IS 
on 10 1 = 1 »Nl 

in mi:=wi*Q ( I ) *x** (i-l> 

GO TO 20 

1 s W 1 - l« 1 ♦ Q (1 ) 

20 continue 

W = W 1 
RPTURN 


007380 
007390 
007400 
007410 
007420 
007430 
007440 
007450 
007460 
007470 
007480 
007490 
007491 
007500 
007510 
Q0752Q 
007530 
007540 
„ 007550 
007560 
007570 
007580 
007590 
007600 
007610 
007620 
007630 
007631 
007640 
007650 
007660 
007670 
„ 007680 
007690 
007700 
007710 
007720 
007730 
007740 
007750 
007760 
007770 
007780 
007790 
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ENH 

00U«LEP«£CJ1>I0N function ZHtM 
IMPLICIT DOUrtLF PRECISION ( A-n f 0-2) 

THTS JUNCTION DEFINES (?) AT DISCRETE X. 

COMM0N/BLK l / e ■» T H * R H 0 * G * M * 8 * W( ) 

COMMON Qd^) 

Ml = NO ♦ 1 
jrsMl^O . 0 

if (x O t:o. o.O) go ro is 
IF (R .Eu # X) GO To 20 
On 10 I=?*Nl 

10 TPpwls' T E k m 1 ♦ (u ( 1 ) *x»* ( i-p) ) * F L 0 A I ( I-I ) 
TFPmi-TFDG] + x# ( -0 .5 ) 

GO TO ?S 
IS TFRM1 a 0(2) 

GO TO 25 
20 CONTINUE 

DO 30 I=2«n1 

30 Tfom] = Term! + (O(I) * X**(I-2)) * FLOAT (1-1) 

25 CONTINUE 

TFPM a -fFpMi 

?PsTERM 

PF TURN 

IFNin 

OOHR-LEPRECISION FUNCTION ZOP(X) 

I mpj I C I I DOUBLE PRECISION (A~M.t.W) 

THIS FUNCTION DEFINES ( /. ) F 0 « 0 I 3 C H E T E X. 

COMMOM/GLK 1 / F * T H » OHO • G 9 p * K * N(i 
COMM O N Q ( i b ) 

N1 a NO ♦ i 

TF PM 1 = 0.0 

TF (X .EG. 0.0) GO TO )s 
IF <R .EG. X) GO TO 20 
00 10 1=3*01 

10 TFP m 1=TEBM1 ♦ ( 6l ( 1 ) ( 1-3) ) *Fi.OA f ( ( J - 1 ) * ( 1 - 2) ) 

TFPMl= TFRMl* (R**<>-x **2 ) »* { -l) . 5 ) 

GO TO 2b 

IS T F O M 1 a 2.0 * U < j) ♦ (1.0 / F ) 

GO TO 2 o 
2 a CONTINUE 

00 30 !=3#Nl 


007800 
007810 
007811 
007820 
007830 
007840 
007850 
007860 
007870 
007860 
0078*0 
007900 
007910 
007920 
007930 
007940 
007950 
0 0 7 9 b 0 
007970 
007980 
007990 
008000 
008010 
008020 
008030 
008040 
008050 
009060 
008070 
008080 
008081 
008090 
008100 
008110 
008120 
008130 
008140 
008150 
U 0 8 1 6 0 
008170 
008160 
008190 
008200 
008210 
008220 
008230 
000240 
008250 
008260 
008270 


ORIGINAL pa™ ro 

OP Pnno J, E IS 

• jor qvalhy 



o n n o 
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C 

C 

C 


c 


c 

c 


c 


c 


3 0 TrpMl^TEKMl ♦ { y ( I) *x** ( t - 3 ) ) * F L 0 A ' ( ( 1 -1 ) * ( 1-?) ) 


008280 

25 CONTINUE 


008290 



008300 

TpPM - -TEPHl 


008310 



008320 

7DP = 1>:9 M 


008330 



0083^0 

RF TURN 


,008350 

End ' 


008360 

DOUPLFHHF-C IS ION f UNCI I ON F'T3 U ,NK,.v) 


008370 

IMPl ICI r DOUBLE PRECISION (A-n.rW) 


008371 

THIS FUNCTION DEFINES THE InIEGhAno FOR OE(M/uQ(M) FOR 

DISCRETE 

Xo 08 38 0 

K=]*2* ( NO + 1 ) AND 1 * ? * ----- • ( NO + i ) . IN THE MATRIX K 

IS THE RO WO 0 8390 

AMD M IS THE COLUMN* 


008400 



008410 

COMmON/BLK I / E 9 T h « RhO • r, , p * kp.^o 


008420 

c 0 M M 0 N 0 (IS) 


008430 



008440 

JF (X .EU. 0*0 .OR. NP .EE. 2) uO TO 0 


008450 

TF| OT] = -FLOAT (( NP-l )* ENP-2 ) ) * ,<**(NP-3) 


008460 

FD TO 10 


008470 

5 TFI OTI = 0,0 


008480 

10 continue 


008490 

IF (X . E 0 . 0.0 .UK, ft|P . L.L . 1) 00 TO ib 


008500 

TFLOT2 = -FI.OAT (f\iP-l ) * X**(Nh-2) 


008510 

on TO 20 


008520 

lb TF| DTP = 0.0 


008530 

20 CONTINUE 


008540 



008550 

\] ? - NO + 3 


008560 



008570 

A 1 = 1,0 + ZP(X)**2 


008560 

PI r ZDP(X) 


003590 

Cl = T F L 0 T 1 


008600 

A-p s 7 DP ( X ) #*2 


008610 

R? s 7 P { * ) 


008620 

r ? = 2 . 5 0 * T F L 0 T 2 


008630 

AT = ( ( 1 . 0+ZP (X) **2) ** (-0 .bO i } 


008640 

E 3 = 7 P ( x ) 


008650 

C 3 - TFLuTP 


008660 

D sr ( 1 . 0 ♦ 7P ( X ) **2 ) **(-3.50) 


008670 

A 1AM s 2.0 * ZP(A) *P/P(X,M) 


008680 

SUM = P7.DP ( X ,M ) 


008690 

A ? A m = 2.0 * 7 0 p ( < ) * p 7 D P ( X * h > 


008700 

Q ! ? A M = P 7 P ( X 9 M ) 


008710 

A 1 A M s -2 p(X) * ( ( i .0 + 7.P ( V) **?) ** (-1 .30) ) * PZP(X.M) 


008720 

p 3 A M = p Z p ( X f M ) 


008730 

DAM = (-7.0) * ((1.0 + ZP(X>**2) ** (-4, SO)} <* 7P(X) * 

PZP ( X »M) 

008740 



008750 

TFPMl = ( (F*TH**3> /l2. 0 ) * 


008760 

1 (D * ((4l*Cl*BlArt) + ( R i 1 * A 1 Am ) - ( A2*C2*82AM) - 

008770 

? (R2*C?*A<AM) ) + (Al*nl*Cl - A2*62*C2> * DAM) 


003780 



008790 


I 






AX-9 



TPpM? = { 0 . 50-“> * IH*Q {M2 ) a*2 ) *' 

008800 



1 (*3«C3*83AM ♦ 8.3*C3*AJAM) 

008810 

c 



008820 



T F (X .E-j. 0 * 0 ) bO To 25 

008830 



TFdm3s (RHO*G) * <PZ ( * * m ) o (- { x** ( np- l ) ) ) ) 

008840 



00 TO 30 

008850 


28 

TfrpM3 = u.O 

008860 

c 



008870 


3 n 

Tprms TERRI ♦TEPM2+TEKMJ 

008880 

c 



008890 



FT3=TFRM 

008900 

c 



008910 



wftiwn 

003920 



FMP 

008930 



nniiOLFPREC T S 1 0 * j FUNCTION FT4(x,nP*m) 

008940 



implicit double precision (a-m,o-z) 

008941 

c- 


This FUNCTION DEF INtSIMtPWANO F OP OF ( K ) /DQ ( rf ) 

FOR DISCRETE X AND 008930 

c 


k=]*2* (AjO+1) AND m=N0+2 

008960 

c 



008970 



COMM()N/tiLK 1 / F ♦ Th , 9t~,0 . G « '-> * 9 *nO 

008980 



COMMON U(13) 

008990 

c 



009000 



M? = NO + 2 

009010 



IF (X .F;j. 0.0 .OH. -VP ,Lt. >) 00 TO IS 

009020 



TF) OT? = -FLUftT(rtP-l) » 

009030 



On TO 2 0 

009040 


\ F 

tfi or? = o.o 

009050 


20 

CONTINUE 

009060 

c 



009070 

c 



009080 



TFPMas (F*t H*0 (N2) ) * ( U .0 + P ( x> *<‘t: )^{-0.b) J ^H( A) 

* TFL0T2) 009090 

c 



009100 



CI4 - T FH ^ 

009110 

c 



009120 



RFTUrtN 

009130 



Fxin 

009140 



OOUPLEPHtCIS I UN r UNCTION FlblX, NR*m) 

009130 



iMPiici f nouALt: precision (A-n,rw> 

009151 

c- 


juy<; FUNCTION DEFINES I M T K GRAND FOP OF t K ) /DO ( M ) 

FOR DISCRETE X AN0009160 

c- 


K-NO ♦ ? AND M= 1*2* "JO+ 1 . 

009170 

c 



009180 



COMMON/ 8 L K ] / E * ) H * HHO . 0 * u * rt 9 < •-; < ) 

009190 



c 0 M M U N w ( 1 •?> ) 

009200 

c 


r • 

009210 



N2 = NO + ? 

009220 



TFQHs( E :* f ( N? ) ) * < ( 0 .3) * { 1 . 0 + /P iX)**2)**I-0*3) 

009230 



1 *{2 # 0*Z a (AJ*P £P ( X * N ) ) ) 

0092*0 

c 



009250 



FT 3 =5 T ER ' 

009260 

c 



009270 



Or TURK? 

009260 



F \i n 

009290 



non 
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C 

C 

c 


c- 


s 

10 

IS 


c 


c 


s 


in 

IS 


c 


OGi IGLEPRLC I SION FUNCTION F I 6 l X « N^ » N ) 
IMPLICIT DOUBLE PRECISION 
THTS FUNCTION DEFINES I N T F (3 R A J 0 FOR 
K=NO+? AND M-=N0 + 2 

COMMON /BLK 1 / E » ! H * R|-*0 * 9 P » K » NO 

COMMON 0 ( i :> ) 

TF»M = (F* I H) * i 1 • U + ZP ( X ) ##£ ) ( u . S) 

FTAsTFWM 

prjURN 

END 

QOHPLFPBtCISlON K UNC X 10 i PZOP(X? w ) 
IMPLICIT DOUBLE precision (A-h#0-/) 
Tmt^ function CALCULATES THE v a l u p of 
T F (X -E<i. 0.0) 00 To 10 
I r (M . GT. 2) GO TO S 
T F o m = 0-0 

GO TO 16 
CONTINUE 
TFPM = -FLOAT ( 

GO TO 1^ 

CONTINUE 
TFDM = 0-0 

JF (M -EC- 3) TERM - -2.0 

CONTINUE 

P7HP * TEPm 

pr thRN 
FMH 

DOURLFPKECISION FUNCTION P*PU*m) 
IMPLICIT DOUBLE PRECISION (A-h.0-7) 
THIS FUNCTION CALCULATES THE VALUE OF 
I f ( X » E 0 . 0.0) GO TO Id 

ir (M -(51 • 1 ) GO 10 5 

T f o m = 0 • 0 

GO to 16 
CO M T I n U E 

TFPM = —FLOAT tM— iJ * /«'•&( ) 

GO TO 16 
tok'T INUE 

-fp-PM = 0.0 

IF ( m . E O * 2) TER** = - 1 » 0 
continue 

p/p - TERM 

Pf TURN 
ENn 

OOtfOLFPPtClSION FUNCTION p2(A*m) 
IMPLICIT DOUBLE PRECISION (A-H,(W) 


009300 

009301 

UF(K)/00(M) FOR DISCRETE X AN0009310 

009320 
009330 
009340 
009350 
009360 
009370 
009380 
009390 
009400 
009410 
009420 
009430 
009431 

0(7) OR P L Dp FOR DISCRETE X, 009440 

009450 
009460 
009470 
009480 
009490 
009500 
009510 
009520 
009530 
009540 
009550 
009560 
009570 
009580 
009590 
009600 
009601 

0 il) OR P7P(X) FOR DISCRE TE009610 

009620 
009630 
009640 
009650 
009660 
009670 
009680 
009690 
009700 
009710 
009720 
009730 
009740 
009750 
009760 
009770 
009771 



o o o 
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C- THIS FUNCTION CALCULATES THE VALUE OF 0(Z)GRPZ(x) FOR DISCRETE X. 
IF (X . EG . 0.0 ) GO TO 5 
TFPM a ) 

GO TO 10 
5 CONTINUE 
TFPM =0*0 

IF (M .EG* 1) TERM = -1.0 
10 CONTINUE 
P7 = TERM 
C 

RETURN ' 

END 

DOURLEPRECISIOH FUNCTION FT 7 < X , I OUM , JDUM) 

IMPLICIT DOUBLE PRECISION <A-H,0-Z) 

FUNCTION TO EVALUATE ARC LENGTH OF DEFORMED BLADDER 

FT7 = OSQRTd.O * ZP(X)*Zp(X)) 
rfturn 

EM n 7 

SUBROUTINE pagehd 

COmmON/LSTaRT/ IRUNNO, I DAT E »NPAGE * UNftME ( 3 ) * T I TLE 1 ( 12 ) » T I TLE2 ( 12) 
DATA NIT, NOT/S, 6/ 

C 

2001 FORMAT {GH1RUN NO. * AH , 42X , 5Hr) A jt *A6/ 

*55X,7HRUN RY f3A6/10?X,A9,12n CLOCK T I ME/ 1 OX , 1 2 A6 ,1 9X , F 1 0 . 3 * 

SEC. CPTIMt/lOX, 1?AH) 

C 

CAM. TIME (DTIME) 

CALL SECOND (CP) 

WP I TE (NOT. 2001) LRUnNO, IDAT t ,UNAME, DTIME, 

* TITLE!, CP* TITLED 

C 

RFTURN 

FMD 

SUBROUTINE MULTw ( A , HZ , NR,A , WRR , NCB , K A * KB Z ) 

IMPLICIT DOUBLE PRECISION (A~H,fW) 

DIMENSION a (KA, 1 ) *BZ (KBZ* 1 ) 

COMMON /LWRK V 1 / W<40) 

c 

C MATRIX MULTIPLICATION.' A * A = Z. 

C USES TWO WORK SPACES. RESULT U) IS PLACED IN R. 

C RZ MUST BE DIMENSIONED LARGE ENOUGH IN MAIN PROGRAM TO CONTAIN THE 

C LAPGFP OF B OP Z. 

c calls forma subroutine zzromb. 

c the maximum size is 

c NPR =40 

C DEVELOPED BY CARL BOOLEY. JANUARY i9E5. 

c last revision by r l wohlen. july ir7?. 
c 

C SURPOUTINE ARGUMENTS 

C A = TNPUT MATRIX, S I ZE ( NR A « NKB ) . 


009780 
009790 
009800 
009810 
009820 
009830 
009840 
009850 
009860 
009870 
009860 
009890 
009900 
009901 
009910 
009920 
009930 
( 0 09940 
009950 
009960 


000100 

000110 

000120 

000130 

000140 

000150 

000160 

000170 

000180 

000190 

000200 

000210 

000220 

000230 

000240 

000250 

000260 



/ 


AX- 12 


c 

87 

s 

INPUT 

matrix 

. S IZE ( NRR * 

c 


s 

OUTPUT 

result 

MATRIX. 

SI 

c 

nra 

= 

INPUT 

NUMBER 

OF ROWS 

OF 

c 

nrb 

a 

INPUT 

NUMBER 

OF ROWS 

OF 

c 

NCR 

= 

INPUT 

NUMBER 

OF COLS 

OF 

c 

ka 

= 

INPUT 

ROW DIMENSION 

OF 

c 

krz 

a 

input 

ROW DI MENS I ON 

OF 


MATRICES A * Z • 

MATRIX COLS OF MATRIX A. MAX=40 • 
MATRICES BfZ. 
a in calling program* 
uz in calling program. 


IF (NP8.GT.40 .OR. NRA.GT.KBZ .OR. NHW.GT .KBZ) GO TO 

DO 40 J=I»NCB 
DO 20 K=WNRB 
2 0 W < K ) = RZ(K*U> 

DO 40 I=i.NRA 
S = 0. 

DO 30 K=1*NRB 
30 S = S + A ( I ,K ) *W (K > 

40 B7(I,J) = S 
RFTURN 

999 C 6 LI 7ZSOMB (6HMULT8 ,NERROh> 

ENtn 

SUBROUTINE ZZ80MB < SU8NAM , NERROR ) 

DIMENSION DFMSSG(B) 

DATA NTT* NOT/5 * 6/ 


NERNOR= 1 
999 


CONTPOI. A COMPUTER RUN AFTER AN ERROR 
IN ANY OF THE FORMA SUBROUTINES. 

ON THF CDC D0Q0 SERIES COMPUTER THIS 
O) PRINT ERROR MESSAGE* INCLUDING 
(?) PRINT ERROR MESSAGE* INCLUDING 
(3> CALL TO NON-EX 1ST ANT 
STOP A NO TRANSFER TO 
COOED RY RL WOrtLEN. SEPTEMBER 1970 
LAST REVISION RY R HRUDa* JflN 1974. 


MESSAGE HAS BEEN ENCOUNTERED 


INVOLVES ... 

SUBNAM AND NERROR* 
SUBNAM AND NERROR* 

rouiine to cause abnormal 
TH t exit card. 


IN 

IN 


PRINTOUT 

DAYFILE. 


SUBROUTINE ARGUMENTS 

SUPNAM = INPUT SUBROUTINE NAME IN WHICH ERROR OCCURRED. 
NEPROR - INPUT ERROR NUMBFR FROM SUBROUTINE WHERE ERROR 


OCCURRED 


300} FORMAT ( 1 H 1 ) 

300? FORMAT ( 1 9H Z280MB - 


ROUTINE ( * A b * 1 1 H ) * NERROh (,I3*1H)) 


WRITE (NOT, 3001 ) 

WPTTE (NOT, 3002) SU0NAM , NERROh 
ENCODE (40,3002,DFNSSG) SUBMam , NERROR 
call remark ( OF MSS (3 ) 

CALL ARNORML 


END 

SUBROUTINE 


8 T ABA ( A Z,B,NRU* NCR, KAZ* KB) 


000270 
000280 
000290 
000300 
0003X0 
000320 
000330 
000340 
000350 
000360 
000370 
000380 
000390 
000400 
000410 
000420 
000430 
000440 
000450 
000460 
000470 
000480 
000490 
000500 
000510 
000520 
000530 
000540 
000550 
000560 
>000570 
000580 
000590 
000600 
000610 
000620 
000630 
000640 
000650 
000660 
000670 
000680 
000690 
000700 
000710 
000720 
000730 
000740 
000750 
000760 
000770 
nnnftl o 



Al-13 



implicit double precision <a-h;o~z> 


000820 


DIMENSION AZ ( K A Z * 1 ) * 8<KB*1) 


000830 

' 

COMMON / LWRKVI / W(40) 


000840 

c 



000850 

c 

TRIPI F MATRIX PRODUCT. H ( TRANSPOSE ) * A * B = Z. 


000860 

c 

A MUST BF SYMMETRIC TO GET CORRECT ANSWER. 


000870 

c 

? WILL BF SYMMETRIC, upper half calculated* reflected to 

LOWER HALF. 

000880 

c 

uses two work spaces, result u) is placed in a. 


000890 

c 

A 7 MUST' BE DIMENSIONED LARGE ENOUGH IN MAIN PROGRAM TO CONTAIN THE 

000900 

c 

LARGER OF A OR Z* 


000910 

c 

CALLS FORMA SUBROUTINE ZZBOMH. 


000920 

c 

THE maximum sizes are 


000930 

c 

NRg = AO 


000940 

c 

NCR = 40 


000950 

c 

DEVELOPED BY w A REnFIELD. MAY 1972. 


000960 

c 

LAST PFVISION BY R a PHILIPPUS. JUNE 1972. 


000970 

c 

' 


000980 

c 

SUBROUTINE arguments 


• 000990 

c 

A 7 = INPUT inner matrix. SIZE (Nk8*NH8) . 


001000 

c 

- output result matrix. S I ZE (NCH »NCB ) . 


001010 

c 

B a INPUT OUTER MATRIX. SIZE (NRpiNCfl) . 


001020 

c 

NPB = INPUT NUMBER OF ROWS OF MATRIX B * SIZE OF MATRIX 

A. MAX=40. 

001030 

c 

NCR = INPUT NUMBER OF COLS OF MATRIX B » SIZE OF MATRIX 

Z. MAX=40. 

001040 

c 

KAZ s TNPUT ROW DIMENSION OF AZ IN CALLING PROGRAM, 


001050 

c 

KR = INPUT ROw DIMENSION OF ti IN CALLING PROGRAM. 


001060 

c 



001070 



nerror=i 

001080 


IF (NRB.GT.4U .OR. NCP.GT.4U .OR. NRR.GT.KA/ .OR. NCB 

•gt.kaz) 

001090 


* GO TO 999 


oonoo 

c 



OOlllO 


DO 20 I=i*NRB 


001120 


no 5 K= 1 * NPB . 


001130 


5 W <K > = AZ ( I *K) 


001140 


on 20 J=1*NCB 


001150 


s = 0.0 


001160 


no io k=i.nrb 


001170 


10 S = S + W(K)*B<K*J) 


001180 


2 0 A 7 ( T * J ) = S 


001190 

c 



001200 


DO 30 J= 1 * NCR 


001210 


DO 25 I=1*J 


001220 


W ( T } = U* 0 


001230 


DO ?5 K=1 ,NR8 


001240 


2S W ( T ) = W(I)4>b(K*I)*AZ(K.J) 


001250 


no 30 i=1*j 


001260 


A 7 { T * J ) = *(I) 


001270 


30 A 7 ( J * I ) = a (I ) 


001280 


RETURN 


001290 

c 



001300 


999 CALL ZZBOMR (oHBl'ABA *N£RRUK) 


001310 


END 


001320 


SURROU T I NE IN'/I (AjZ*N.KP) 


001330 



o n 
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C 

C 

C 

C 

C 

C 

C 

C 

C 

c 

c 

c 

c 

c 

c 

c 

c 

c 


c 

c 

c 


IMPLICIT DOUBLE PRECISION (A-h,0 -Z) 

OTVENS ION A ( 1 ) * 1 ( 1 ) 

COMMON /LwPKVi/ 0(20) » 0ETH(20) 

COMMON /L WRKVH/ IX (20), 8(20) 

DATA N IT * NOT AS *6/ 

MATRIX INVERSION (A**-l = Z). BORDERING METHOD* 

THE DETERMINANT PATIO DET(I+1) / OET(I) IS PRINTED. DET ( I ) IS THE 
DETERMINANT OF THE FIRST I BY I SUB-MATRIX OF A. 

THE INVERSION CHECK Z * A IS CALCULATED AND PRINTED* 

matrtcfs a*z may shape same core locations. a check is invalid). 

CALLS FORMA SUBROUTINES P AGEHO * Z/ROMB . 

The maximum size is 

N s 20 

DFVELOPED BY BOB DILLON. FEBRUARY 1968. 

LAST REVISION BY J ERNST, OCT 1973. 

SUBROUTINE ARGUMENTS 

A = tnjPUT matrix to be inverted. S I Z E { N * N ) , 

Z = OUTPUT RESULT MATRIX, $J2E(N,N). * 

N = INPUT SIZE OF MATRICES A*Z. MA X=2U. 

KP = INPUT ROW DIMENSION OF A , Z IN CALLING PROoHAM, 


2000 FORMAT 

2001 FORMAT 

2 on? format 
<& 

2003 FORMAT 
•a 


(// lOX*lU(7X*lH(*l2fiH))) 

(// 1 0 X • 4 S H S Up P OU T I NE INV1 HAS CALCULATED THE DATA BELOW 
///10X,44 hTHE DETERMINANT PATIOS DET(I+l) / DET ( I ) ARE 
// ( I3X* 10011.3) ) 

( /// 1 OX , 3 ?H 1 HE INVERSION CHECK GIVES 

///10X,2SHTHF DIAGONAL ELEMENTS are // ( 13X.8014.6) ) 

{// 10X.3SHTHF MAXIMUM OFF-DIAGONAL ELEMENT IS 
Dll. 3* 2A* 4 H a I < 13, 1H, lb, 1H) ) 


IF (N .GT. 20) GO TO 999 


NERROR=l 


DO 160 I=2,M ' 

160 t x ( r ) = r 

INVERT FIRST NON-ZERO ELEMENT I N FIRST COLUMN. 

DO 190 I = 1 * N 

IF ( A { I ) . NE • 0.) GO TO 220 

190 .continue 

' NERROR=2 

GO TO 999 


START INVERSION WITH ROW I. 
220 OFT» U ) = A < I > 

Z ( T ) = 1 • / A ( I > 

IF (N ,E0. 1) return 
C 

IX { I ) = i 

I X { I ) = I 

C BQPDFRING loop. 


001340 
001350 
001360 
001370 
001380 
001390 
001400 
001410 
001420 
001430 
001440 
001450 
001460 
001470 
001480 
001490 
001500 
001510 
001520 
001530 
001540 
001550 
001560 
001570 
001580 
001590 
001600 
001610 • 
001620 
001630 
001640 
001650 
001660 
001670 
001680 
001690 
001700 
001710 
001720 
001730 
001740 
001750 
001760 
001770 
001780 
001790 
001800 
001810 
001820 
001830 
001840 
001850 
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00 630 L = 2*N 
K = L 

LI = L - 1 
250 S = 0. 

M I XL = KR * (1A(L) 
LL * IX(L) + MIXL. 

DO 450 I = 1 f L 1 

M I X I r KP * ( IX ( I ) 

LI = I X < L ) ♦ M 1 X I 

B<I> = 0. 

cm =0. 

no 440 J=l*Ll 

MIXJ = KR * (IX (j) 
IJ - I X C I ) + Mlxj 
JL = IX(J> + MIXL 
R { T ) = B ( I ) - 7(1 J)* A { JL) 

JI = IX(J) + M I XI 
LJ = I X ( L ) + MlXJ 
440 G ( T > = G ( I ) - A(LJ)* Z-(JI) 

450 S x S + A <LI ) * 8 ( I > 

AL = A (LL) + 5 

IF (A (LL) .EG* O.) GO TO 460 
Aj PAR x DABS (AL / A (LL) ) 

GO TO 490 

480 ALRAR = OABS (AL) 

490 IF (ALBAR .GE* .10-6) GO TO :>S0 

INTERCHANGE ROWS AND COLUMNS. 

K = K ♦ 1 

IF (K .GT. N) C 7 0 TO S4 0 
I* L = I X (L ) 

T X f L > = I X { K ) 

IX (K) = IX L 
GO TO 2b0 

540 IF (ALBAR . GE . .10-8) GO TO 550 

GO TO 99R 

550 7 (1 I. ) = 1 • / AL 

l)FTR(L) = AL 
00 570 I=l*Ll 

IL = IX(I) ♦ Ml XL 
LI = IX (L) + NR * ( IX ( I ) - i > 

7 ( TL ) a BU) * 7 (LL) 

7 (I. T ) = G ( t) * l (LL) 

00 570 J - 1 § L 1 

u s rxm + kp * (ix(j) - i) 

570 7 ( T J) - 7 ( I J) + G(J) * 7 (I L ) 

630 CONTINUE 

COMPUTE INVERSION CHECK 2* a. 

X OFF = 0.0 


001860 

001870 

001880 

001890 

- 1) 001900 

001910 

001920 

- 1) 001930 

001940 

001950 

001960 

001970 

- 1) 001980 

001990 

002000 

002010 

002020 

002030 

002040 

002050 

002060 

002070 

002080 

002090 

002100 

002110 

002120 

002130 

002140 

002150 

002160 

002170 

002180 

002190 

002200 

N£RR0R=3 002210 

002220 
002230 
002240 
002250 
002260 
002270 
002280 
002290 
002300 
002310 
002320 
002330 
002340 
002350 
002360 
002370 
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nn 720 1=1 , n 

DO 710 J= 1 , M 
X = 0.0 

K Jft = KR * ( J-l ) 

DO 703 K=1*N 

s I + KR* (K- 1 ) 

KJ = K + KJA 

703 X » X * Z ( IK) * A (KJ) 

IF (I .NE. J) GO TO 70b 
0(T> = X 

GO TO 710 

70S IF (DfiBS(X) .LT. UABSUOFF)) GO TO 710 
XOFF * X 
TOFF = I 
J OFF = J 
710 CONTINUE 
7?0 CONTINUE 

PRINT THF DETERMINANT RATIO AND INVERSION CHECK • 

CA( {. PAGEhi) 

W9TTF (NOT, 20 00) (JC, JC=1«10> 

WRTTE (NOT, 2001 ) (UETR(I), I=1,N) 

WRITE (NOT, 2002) ( G ( I ) * I = 1 , N ) 

WPTTF (NOT, 2003) AQFF, IOFF? JOFF 
RETURN 

999 CALL 7ZBOMB (6HINV1 »NERR0H) 

FNO 

SUBROUTINE MULT ( A , 0 , / , NRA ♦ NRH • NCR , KR A *KRB ) 

IMPLICIT DOUBLE PRECISION (A-H,(W) 

DIMENSION A(KRA*l)f B ( KRB » i ) , Z(KRA»1> 

C 

C -MATRIX MULTIPLICATION. a * B = 7. 

C DFVFl OPED BY PL WOrtLEN « FEBRUARY 1966. 

C LAST REVISION BY R L wOHLEN • JuLY 1972. 

C 

C SUBROUTINE ARGUMENTS 

C A = TNPUT MATRIX. SIZE (NRa »N«B) « 

C ft a INPUT MATRIX. SIZE (NRB, NCR) . 

C l OUTPUT RESULT MATRIX. SIZE (NRAtNCB) . 

c NR A a INPUT NUMBER OF ROWS OF MATRICES A ,Z. 

c Z OUTPUT RESULT MATRIX. S IZt (NRA ,NCB ) . 

C NRB - INPUT NUMBER OF ROWS OF MATRIX H* COLS OF MATRIX A. 

C NCB “ INPUT NUMBER OF COLS OF MATRICES B*Z. 

C KP A a INPUT ROW DIMENSION OF A*Z IN CALLING PROGRAM. 

C KRB a TNPUT ROW DIMENSION OF B IN CALLING PROGRAM. 

C 

DO 20 1 = 1 * NRA 
DO 20 J = 1 * NCR 
S a 0. 

DO 10 K=1,NRB 
in S a S + A ( I 9 K) *6 (K * J) 

20 7 ( T * J) = S 


002380 

002390 

002400 

002410 

002420 

002430 

002440 

002450 

002460 

002470 

002480 

002490 

002500 

002510 

002520 

002530 

002540 

002550 

002560 

002570 

002580 

002590 

002600 

002610 

002620 

002630 

002640 

002650 

002660 

002670 

002680 

002690 

002700 

002710 

002720 

002730 

002740 

002750 

002760 

002770 

002780 

002770 

002790 

002800 

002810 

002820 

002830 

002840 

002850 

002860 

002870 

002880 

002890 
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C 

C 

C 

C 

C 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 


c 


c 


c 

c 

c 

c 

c 

c 


rfturn 

END 

SUBROUTINE MlJLTA ( AZ *NRA *NHH »NC8*K AZt KB > 

' IMPLICIT DOUBLE PRECISION <A-H,Q-*/> 

DIMENSION AZ(KAZ*I>* B ( KB * 1 ) 

COMMON / LWRKV1 / WUO) 

matrix multiplication, a * e = z. 

USES TWO WORK SPACES. RESULT (Z> IS PLACED IN A • 

A 7 MiJCT RE DIMENSIONED LARGE ENOUGH IN MAIN PROGRAM TO CONTAIN THE 
L6RGFR OF A OR Z. 

CALLS FORMA SUBROUTINE /ZBOmh. 

THE MAXIMUM SIZE IS 

NP» s SCO 

DEVELOPED BY c S BOOLEY. JANUARY 196b. 

LAST REVISION BY R F HRUDA • , JUNE 1972. 

SUBROUTINE ARGUMENTS 
A 7 - INPUT MATRIX. SIZE (NRA *NKB) , 

= OUTPUT RESULT MATRIX. S I Zt < NR A , NCB ) . 

B = INPUT MATRIX. SIZE <NRB«nCR) 

NR A = INPUT NUMBER OF ROWS OF MATRICES A»Z, 

NRB = INPUT NUMBER OF ROWS OF MATRIX Rf COLS OF MATRIX A. MAX=500. 

NCB = INPUT NUMBER OF COLS OF MATRICES 6*2. 

KA? = INPUT ROW DIMENSION OF AZ CALLING PROGRAM. 

KR - INPUT ROW DIMENSION OF B IN CALLING PROGRAM. 

NEHKORel 

IF (NPB .GT. BOO) GO TO v99 

DO 40 I - 1 » NR A 
DO 20 Ksl.NRH 
?0 W (K) = AZ ( I ,K ) 

DO 4 0 J = 1 * N C B 
S - 0.0 
DO 30 K=1 . MHd 
30 S ~ S + w(K)*ti(KiJ) 

40 a 7 ( I * J > - S 
RF TURN 

9 99 CAM. ZZBO w B (6HMULTA *NERROR) 

END 

SLJOPOUTINE WRITE { A * NR * NC * I N a mE * K R ) 

ImpiICIT DOUBLE PRECISION (A-H.O-Z) 

DTMFNSIDN A ( A R » 1 ) 

DATA N IT*N0T/S*6/ 

WrITF MATRIX QF WE A{. NUMBERS ON PAPER. 

RF <)U T °FS 123 COLUMN (MInIM'JM) PRINTER. 

UP TO 10 DA T A FIELOp PER LINE. PRINTS ONLY NON-ZERO FIELD ROWS. 
CALLS FORMA SUBROUTINE PAGE HI). 

COOFD BY RL WOHLEN. DECEMBER 1968. 


002900 

002910 

002920 

002930 

002940 

002950 

002960 

002970 

002980 

002990 

003000 

003010 

003020 

003030 

003040 

003050 

003060 

003070 

003060 

003090 

003100 

003110 

003120 

003130 

003140 

003150 

003160 

003170 

003160 

003190 

003200 

003210 

003220 

003230 

003240 

003250 

003260 

003270 

003260 

003290 

003300 

003310 

003320 

003330 

003340 

003350 

003360 

003370 

003380 

003390 

003400 

003410 
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C LAST RPVISION 0Y R HRUOa, .MOV 1*73. 

C 

C SUBROUTINE ARGUMENTS (ALL INPUT) 

C A = MATRIX TO BE PRINTED. SIZE (NR. NO. 

C NR = number of rows in MATRIX A. 

C NC = NUMBER OF COLS IN MATRIX a. 

C IN AMR s MATRIX IDENTIFICATION. lA6 FORMAT). 

c kr = row dimension of a in calling program, 
c 

2010 FORMAT <//15H OUTPUT MATRIX A6.2X 1H(I4,2H X 14, 2H ' ) // 
* 10X,10(7X,1H( I 2 , 1 H ) ) / ) 

2020 FORMAT (//ISH OUTPUT MATRIX A6,2X 1H(I4»2H X l4,2H ) 

'* 3X, 9HCONTINUED /✓ 1 OX • 1 0 ( 7X , 1 H < I2*lrt))/) 

2030 FORMAT (1X.2I5.2X, 10O11.3) 

2040 FORMAT ( 14H0END OF WRITE.) 


C 


c 

c 


c 


PULL UP A NEW RAGE FOR MATRIX ANf) PRINT MATRIX NAME. 
CALL PAGEHD 

WRITE (NOT, 2010) INANE .NR.NC » (L*L=1 , 10) 

NLINE = 0 


DR 50 1=1 »NR 
N7FR0 = 0 
JS = 1 

10 JF = JS+O 

IF (JF .GT. NC) JE=NC 
SFF IF ELEMENTS ARE 2£R(). 
no ?0 J=JS,JE 

IF ( A ( I , J ) « NE . 0.) GO TO 30 

20 CONTINUE 
GO TO 40 

3 0 NLTNE = NL.JNE+1 

IF (NLINE .LE. 4*0 GO TO 3b 
CALL PAGEHD 

write { NO T ,2020) I NAME * NR * NC , ( L .L= 1 , 1 0 ) 

NLTNE = 1 

35 WRITE (NOT, 2030) l . JS , ( A ( I * J ) , J = JS.JF) 

N7FR0 = 1 

40 IF (JE .EO. NC) GO T(j 50 
JS = JS+10 
GO TO 10 

SKIP A SPACE BETWEEN EACH ROW IF THERE APE MORE THAN 10 COLUMNS 
AND SOMETHING HAS BEEN WRITTEN. 

50 IF (NC.LE.10 .OR. NZERQ.EO.O .OR. I. EG. NR) GO TO 60 
NLTNE = NLINE+1 
write (NOT, 2030) 

60 CONTINUE 


WRTTE (NOT, 2040) 

RFTURN 

END 



003420 

003430 

003440 

003450 

003460 

003470 

003480 

003490 

003500 

003510 

003520 

003530 

003540 

003550 

003560 

003570 

003580 

003590 

003600 

003610 

003620 

003630 

003640 

003650 

003660 

003670 

003680 

003690 

003700 

003710 

003720 

003730 

003740 

003750 

003760 

003770 

003780 

003790 

003800 

003810 

003820 

003830 

003840 

003850 

003860 

003870 

003880 

003890 

003900 

003910 

003920 
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SUBROUTINE ZERO ( 2 * NR « NC * KR ) 

IMPLICIT DOUBLE PRECISION (A-h,0-Z> 

0 1 MANSION Z(KR*i) 

c 

c generate a matrix of zeroes. 
c coded by RL wohlen. feb 1965. 
c 

c SUBROUTINE arguments 

C Z r OUTPUT MATRIX GENERATED. S1ZE(NR*NC). 

C NP = tnput number of rows in matrix z. 

C NC = INPUT NUMBER OF COLS in MATRIX Z. 

C KR r input row DIMENSION of matrix z in calling program. 

c 

DO 10 1=1* NR 
00 10 J= 1 ♦ NC 
10 Z(T*J) = 0.0 
PFTURN 
END 

SUBROUT IMF START 

DI M ENSIQN mONThn (12) .MONTHL (12) 

COMMON /LST ART / IRUNNO * ID A TE * MPAGE t UN A ME < 3 ) • TITLE! < 12) »TI TLE2 02) 
DATA N I T « NO r/6 * 6/ 

DATA MONTHN/2H01 *2H02»?H03*2o04 t2HO5*2H0fe, 

* 2H0 7 *2H0S*2HO9*2H10 *2Hl 1 *2H12/ ♦ 

* M0NTHL/2HJA »2HFE *2HMRt2HAP*2HMY»2HJN* 

* 2H JL * 2HAU » 2HSE * 2 HOC * ?HN0 * 2HDE/ 

C 

C READS INPUT CARO 1 FOR IRUNNO* UNAME. 

C CHECKS IRUNNO FOR STOP (I.E. IF IRUNNO * STOP* PROGRAM WILL 
C BF STOPPED). YOU SHOULD HAVE A STOP CARD (THE WORD STOP PUNCHED 
C STARTING IN COLUMN i) AFTER YOUk REGULAR DATA DECK. 

C IF IRUNNO lb NOT EQUAL TO STOP* THE SUBROUTINE CONTINUES AS FOLLOWS. 
C READS INPUT CARD 2 FOR TI TLEI . 

C READS INPUT CAKU 3 FOR TITLE?. 

C sfts WPAGE - 0. 

C INTERROGATES COMPUTER TO DEFINE DATE AS AN A6. 

C TMTFROGATES MACHINE FOR THE Tl^t OF DAY AND THE CPTIME 

c and pptnts these items on a shelt of Tnt output every 
c time this routine is CALLFID. 
c 


c 

INPUT 

ORDER 


c 

I PUNNO * UNAME 

FORMAT ( A6 ♦ 4 X 3Ab) 

c 

tttlfi 


FORMAT ( 1 2 A 6 ) 

c 

titl^? 


FORMAT (1?A6) 

L 

c 

PFF INI 

TIONS 


c 

IRUNNO = 

RUN NUMBER. ( A6 FORMAT) 

c 

IOATE = 

DATE. 

(Ah F ORMAT ) . 

c 

NP AGE = 

PAGE NUM0LR. 

c 

UNAMF = 

USERS 

NAME. ( 3 a6 FORMAT) 

c 

TITLE! - 

FIRST 

TITLE. ( 1 ? A 6 FORMAT) 

c 

tttle? = 

SECOND 

TITLE. (12A6 FOhmaT) 

c 

c 

dedicated 

TO G. 

MOROSOW. 


003930 
003940 
003950 
003960 
003970 
003980 
003990 
004000 
004010 
004020 
004030 
0 04 040 
004050 
004060 
004070 
004060 
004090 
004100 
004110 
004120 
004130 
004140 
004150 
004160 
004170 
004180 
004190 
004200 
004210 
004220 
004230 
004240 
004250 
004260 
004270 
004280 
004290 
004300 
004310 
004320 
004330 
004340 
004350 
004360 
004370 
004380 
004390 
004400 
004410 
004420 
004430 
004440 
004450 
004460 
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C CALLS COC-6000-SERILS COMPUTER STANDARD ROUTINE* DATE , SECOND, TIME. 004470 


C CALLS MARTIN MARIETTA SPECIAL ROUTINE PPTIM. 004480 

C CODED FOR YOUP CONVENIENCE BY YOUR FRIENDLY METHODS GROUP. DEC 1968.004490 
C LAST REVISION BY H HRUOA, MAY 1974. 004500 

C , 004510 

1001 FORMAT ( A6 , 4X 3A6) 004520 

100? FORMAT (I2a6> 004530 

200? FORMAT < 1 H 1 fa(/> 55X 1 OHT I ME SHEET / 3bX 45(1 h-) // 004540 

* 38X 30HCUHRENT TIME OF DAY IN H,M,$ = AlO // 004550 

* J8X 26HT0TAL CPTIME USED TO NOW = F1U.3, 9H SECONDS. / 004560 

* 38X 2 fa H T 0 T A L PPTIMfc USED TO NOW = 16, 4X, 9H SECONDS. ) 004570 

2003 FORMAT (36H1END OF INPUT DATA HAS BEEN REACHED.) 004580 

5001 FORMAT ( 3 ( 1 X » A2 ) ) 004590 

500? FOPMAT ( 3 a 2 ) 004600 

C 004610 

CALL TIME (DTIME) 004620 

CAJ.L SECUNO (CTXME) 004630 

WRTTE (NUT, 2002) DT IME , CT I ME , I PT I ME 004650 

GAEL DATE ( IUATE ) 004660 

DECODE (9, 50 01, I DATE) 1*1, ID, I Y 004670 

DO 20 1=1,12 004680 

IF ( IM.EU.mONTHN ( I ) ) GO TO 30 004690 

?0 CONTINUE 004700 

30 IM = MONThl ( I ) 004710. 

ENCODE (fa, 500 ?, IOaTE) ID,Im,iy 004720 

C 004730 

RF AO (NIT, 1001) I RUNNO , UN A ME 004740 

IF (IRUNNO .NE. 4HSTOP) GO TO 10 004750 

WRITE (NOT, 2003) 004760 

STOP 004770 

C 004780 

10 RF AD (NIT,10U?) IITLFI 004790 

RF AD (NIT, 1002) TITLF.2 004800 

NPAGE = 0 004810 

RETURN 004820 

END 004830 


APPENDIX - A2 
DYNAMIC ANALYSIS PROGRAM 


A2-1 


Si Hi ShPOOOOQ*F2«MA IN 

1 c 

2 C OVERLAY MAIN PROGRAM TO CALCULATE SHUTTLE EXTERNAL TANK SLOSH MOD^S* 

3 C DEVELOPED a Y to fiENFIELD, c BODLEY,'R PHiLjPPUS, R WGHLEN , JULY 19 ? 3, 

9 C LAST REvISJOn BY R A PH i L I PPUS • JULY 19 7 9# 


6 C •*••*#***»••*•«,«**«*•«««•**«*«» 

7 c input da t * Read in this program. 


a c io call start 

9 C CALL COtiENT 


13 

C 

IFINI T * TAPElu 

• 

1 1 

C 

15 I OP T 



12 

c 

IF 

( iOPf 

.to. 

AHSTART ) 

1 3 

c 

IF 

CJOFT 

• eq. 

6HGNXYZ ) 

19 

c 

I F 

( I OP T 

• eq . 

6HFINEL ) 

15 

c 

If 

l I OP T 

.EQ. 

6HXTRAMK) 

16 

c 

IF 

{ I OP T 

• E Q . 

6HMQ0ES > 

17 

c 

♦ IF 

C 1 0 P T 

* E Q , 

6HM0DES ) 

18 

c 

IF 

( l 0 P T 

.to 4 

6HREDUCE ) 

1 9 

c 

If 

( I 0 P T 

• EQ. 

ghsavemk ) 

20 

c 

IF 

( lOPT 

• E U * 

6MSUBSTR ) 

2 1 

c 

IF 

< I OP T 

.tQ, 

6HPL0T ) 

22 

c 

GO 

TO lb 



23 

c 






FORMAT ( 2 A6 | 
FORMAT ( A 6 1 

GO TO iO 

CALL GN*Y 22 ISLE SU8KT FOR InPUT) 

call finel i I nput data from gnxyzzi 
call PXTRa 4SeE SUBKT FOR INPUT) 
CALL GNIMd t S t E SUBKT FOR INPUT) 
CALL OYMODE (SfcE SUBRT FOR IwPUT) 
CALL REDUCE (SEE 5UBHT FOR INPUT) 
CALL SAVE (See SUB R T FOR INPUT) 

CALL SUbSTR (SeE SUBRT FOR InPUT) 
CALL SPLTI (See 5UBRT FOR INPUT) 


29 

2b 

26 

27 

28 

29 

30 

31 

32 

33 
39 
3b 

36 

37 
36 
39 
90 
9 i 

92 

93 
99 
9b 

96 

97 

98 

99 
80 

51 

52 

53 
59 
5 5 

56 


COMMON / RWTaPS / NUTfcL *NUTXY2 »NUTLT *NUTST *NDTNX#NUTNX tNUlBX 

COMMON / bTapeF / NUTMF , N U T K F , N U T if * N U T 2 F » N U T 3 F 

COMMON / BTapeT / NUThjT ,NUT* r ,HUT i T ,NUT2T jNUTjT ,NUT9l ,NwI5T 

COMMON / btAPEM / NUTMM ,NUTKM »NUTTM ,NUTPM ,NUTfM »NUT 1 M , NU T 2 M , NU T j N , 

* NUT9M ,NUT 511 ,NUT6M ,NUT7M 

COMMON / BTAPEU / NUT^D t NUTLO»NUTDO,NUTlD,NUT20»NUT3D»NoT9D 

COMMON / bTAPEB / NUTKB l NUTa6#NUTPa # NUTlb $ NUT2B,NUT36>NuTHB. 

* NUT5B ,NUT6d |NUT7B 

common / btapec / nuthc ,nut*c ,nuttc ,nut ic ,nut2C ,nuT3c , nuthc .nuTuc s 

P N U T 6 C 

COMMON / BTAPER / NUTMR,NUTNK*NUT?8*NUTifl t NUT2K|NUT3R,NUl9R.NUluK, 

* NUT6H 

COMMON / 6 T ApE S / NOTMS,NUTaS*NUTTS 

COMMON / bTAPEP / NUTMP .NUTKP ,nOT TP ,NUt 1 P ,NUT2F »NUT3P ,N0T9P jNUf^P 
COMMON / BTAPEA / NUTPA » NUTF A »NUT I A »NUi 2A |NUT3A »NUT9 a ♦ N uT bA , 

* NUT6A,NUT7A 

COMMON / RESiaP / N85VTI 

COMMON , RTRamS / IFTRAN 

c 

data ni r , not/s » 6/ 

c define read, write tape units for all overlays, 
data NUTELiNuTXYZ / 

* 2 9 , 30 / 

data NUFLT »NUTST ,NUTrtX»NUTKX,NU]BX / 

* 31 , I » 2 , 26, 27/ 

c DEFINE BUFFER In, OUT TAPE UNITS FOR FINELO. 

Data NUiMF,NgTKf : »NUTlF,NUT2F,NUT3F/ 

* 21 * 22 , ( [ , 12 , 13 / 

c DEFINE buffer IN, OUT TAPE UNITS FOR EXTRA M , K OVERLAY. 

DATA NUTMT » N u T K T ,NUT it ,NUT2T , N U f 3 T , N u T 4 T iNUTbT / 

! 2i » Z2 > U , 12, 1 3 » 19, lb / 

C DEFINE BUFFER In, OUT TAPE UNITS FUR YMOOE2. 

DATA NUIMM,NUTKM,NUTpM,NUTTM,NOiFM»NUTtM»NOT2M»NUT3M,NUT9h ,NUT&p t / 
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57 

5b 

59 

60 
61 
*2 
63 
6** 

65 

66 
67 
6 B 

69 

70 

7t 

72 

73 

74 
7b 

76 

77 
76 

79 
60 
8 i 
82 

83 

84 

85 

86 
87 
66 

89 

90 

91 

92 

93 

94 

9b 

96 


' 21 » 

data nuj6m/ 

• 16/ 

DATA N 0 T 7 M / 


10 ? 

lOb 

109 

no 
1 1 1 
112 
1 1 3 


22 , 23, 2% 25, 


1 1 


12 


13, l 4 * 


!&/ 


* 17 / 

DEFINE buffer In, OUT TAPE UNITS FOR STATIC DEFLECTION CALCULATIONS, 
Data NUTKO ,NUTUD ,NUTDD ,NUT ID ,NUT2D »NdT3D ,NUT4D / 

* 22, 2 1 f 2J, 11, 12, 1 3 , 1*4/ 

°lfine suffer in»uut tape units for buckling load calculations* 

Data NUTKB»NUT8B»NUTPb,NUT 18,NUT28»NU73B,NUT40»NUT5B ,NUT68 ,NUT7 a 
22 t 2i» 23, 11, 12, J3, J4i 15, 16, 

U ¥ IT r f hi . flliT T A O r 1 1 V T C C- n W rfl luC T - . -r IHW .,ur A 1 1 Wt ^ r r- a* / 


I NG LOAD CALCULATIONS* 
38 ,NUT4B ,NUT5B ,N0T68,NUT7 a / 

* 22 1 2i| 23, 11, 1 2 » 13* 14, 15, 16, 1 7 / 

DEFINE BUFFER In # 0Ut Tape Units for CONSTANT VOLUME FLUID EL E M 

CALCULATIONS* 

Data NUTMC f NUIKC ,NUTTC ,NUT1C ,N0T2C ,NDT3C ,NUT4C ,NUTSC ,NUT6C/ 

* 21, 22, 24, 11, 12, 13, 14, 15, 16/ 

DEFINE BUFFER In, OUT TAPE units for REDUCING CALCULATIONS, 

Data NUTMR,NUTKR»NUTTR,NUTlR,NUr2R,NUT3R,NUT4R»NUT5R,NUT6R/ 

21 , 22, 24, 11* 12, 1 3 * 14, IS, 16/ 

DEFINE 

data 

22 , 24 , II, 1 2 1 

In, OUT 


r 1 f M v • 9 " • * r\ | nv « * |itv f in | nv 1 jrr | nv | on/ 

21, 22 1 24, 11* 12, 1 3 * 14, IS, 16/ 

JUMPER 1N,0UT TAPE UNITS FOR SUBSTRUCTURE CALCULATIONS* 
NUrMP,NUTKP l NUTTP,NUTlP,NUT2P,NUT3P,NUT4P,NUT5P/ 

*. . _ 21, 22, 24, 11, 12, 13, 14, 15/ 

OEMNE BUFFER In, OUT TAPE UNITS FOR PLOTTING* 

data nutpa,nutfa,nutia»nut2a,nutja*nuT4A»n4 

1 23 , 24 , u, 12, '13, 14* 

. jUFFER In, OUT TAPE UNITS FOR SAVEMfc* 

DATA NUTMS INUTKS , NUTTS/ 

21. 22. ? a / 


,uur TAPE UNITS FUK PLOTTlNvi* 

rFA,NUTlA,NOT2A,NUr3A»NUT4A,NUT5A,NUT6A,NOT7A/ 

* 2 3 , 24 , U, 12, '13* 14, 15, 16, 17/ 

- define buffer in, out ~ ■ — — 

DATA NUTMStNUTKS, 

* 21, 22, 

; define reserve ta^ls? 

Data NRSVTI 

1001 FORMAT { 1 2 A 6 ) 
toio FORMAT l 10* I 


24/ 
/ 28 


10 


15 


IF (NRSVTI 

call start 
call coment 
Read 


10* 1 5 1 

* GT* 01 RE 6' I NO NRSVTI 


call cqment 

Read <NIT,lCOtl I F X N I T ,TAPE 10 
IF. UFI^IT *E0* 6HINIT1L) CAL 
IFTRAM ■ 0 


I f x n i t , tape 10 

6HINIT1L) CALL INTAPE 1 NRSVT I , T APE I D ) 


97 

If 

( I OPT 

•E0.6HSTART ) 

GO 

TO 

10 

96 

if 

( I OP T 

•EQ.6HGNXY2 ) 

GO 

TO 

IS 

99 

If 

UOPT 

• £ Q * 6HF I NEl » 

GO 

TO 

20 

100 

IF 

( I OP T 

• E0.6HXTRAMO 

GO 

TO 

30 

101 

IF 

< jOPT 

• E6},6HM0DES ) 

GO 

TO 

40 

102 

If 

< I OPT 

•E0*6HM£CHE0 J 

GO 

TO 

70 

i 03 

if 

( I OPT 

« E^ « 6HR£DUC£ ) 

GO 

TO 

80 

l 04 

IF 

( I OP T 

, EG • 6HS A VEMK ) 

GO 

TO 

9 0 

1 05 

IF 

( I OP T 

• £<ij * 6HSUB5TR ) 

GO 

TO 

100 

1 06 

if 

I I OPT 

* EG « 6HPL0T > 

GO 

TO 

Ho 


origin al page is 
0F fOOR QUAlin 


la 

20 

30 


00 TO 999 
CALL UNXY22 
00 TO 15 
CALL OFINEL 
00 TO 15 
CALL PXTRA 


NtRRURa 1 



1 IH 
\ lb 

I l 6 
1 1 / 

i i a 

I I v 
122 
121 
122 
123 
i 2 H 

125 

126 
127 
126 

129 

130 


GO TO *S 

ho call gnimd 
call otmode 
GO TO IS 
70 CALL MECHEQ 
GO TO lb 
00 CALL reduce 
GO TO 15 
?D CALL SAVE 
GO TO *S 
100 CALL SUdSTR 
GO TO IS 
| 1 0 call OPlTIO 
go TO is 

c 

?<?? call zzoomb iahpsuosh *nerrqr ) 
end 


i s* P k T Pi *6 ASIC 
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siRiSHPQoaoa*Fi ♦basic 


1 

2 

3 

H 

5 

6 
7 


SUBROUTINE BaSIC (XYZ »JUOF , £ UL • K RX , X C * , K R J , KC J * * R£ • A C fc , 

* NU TEL • NU TX YZ J 

C 

C SUBROUTINE TO REA^ BASIC FINEL DATA FROM caKO INPUT AND WRITE ON 
C NUTEU AND' im U T X Y / « 

C DEVELOPED BY W a eENFlELUf APRIL 1974. 


u 
9 
*C 
1 1 
l 2 

n 


DIMENSION XYZ(K*X,l), JOOF UKJ,1 | , EULtKKE t l ) 
DIMENSION IOaTA <141 

c 

data nit»not / s,6 / 
iogi format u3A6,a2) 

REMIND NUTtL 


I s 
is 
i 6 
i 7 
ia 
IV 
20 
21 
22 
23 
lH 
z s 
26 


C READ JOINT XYZ COORDINATE MATRIX. 

10 CALL REaO (xT2 8 NJ ? nCX»ARX # .XCXI 
C READ JOINT DEGREE OF FREEDOM MATRIX, 

Call k£adjm ( jdof *nrj ,ncj#krj,kCji 
C READ JOINT EULEk ANGLES, 

CALL READ (£UL t NRE t jvCE t KRfc,KCL) 
c READ ONE Card WITH NAMEl FOR FInEL. 

C READ DATA CAROS FOR fXlAL OR BAR OR T^NGL, ETC. 
15 READ <Nl T , lOOl ) idata 
WRITE 'CnUTEL, l 00 I I IDATA 
IF t i 0 A T A ( l),E0.6H RETURN) GO TO 100 
GO TO IS 
C 


2? 

2b 

29 

30 

31 

32 

33 
3 4 


ICO REWIND nutel 

REWIND NUTXV2 

WRITE IwUTXYi) NJ,NCX ,NRJ»NCJ t NRE,NCE 
WRITE <NOTXrz> l ( JOOF £ I » J ) , I <* 1 ,NRj) t J» | , NC J ) 
WRITE INUTXY2) U XYZU.J>»1*1» n J 1 » J» l , NCX ) 
WRITE ( NOT X Y 2 J U E 0 L l I > J) » l * l » N K E I » J = | , N C E > 

rewind notxyz 

Return 


3 S C 

36 £nd 


&PRT Ft.GNlMQ 
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S 4 R I SHPOOOOQ*F | .ON I HO 


subroutine gniho 


INPUT DATA read in THIS PROGRAM. 

CALL TREAD (INITIAL 0 1 SPL MuOfcS} 
RETURN 


8 
9 
I 0 
I I 

u 

I i 
l H 


COMMON /aTAPEM^ NUTM,JjIUTK,NUTTR*NUTZ»NUTI,NUT2#NUT3»nUT9»NUTS, 
* NUT6,NUT7 

dimension vcizoooi* lvuzoooi 
data KV / I 2 c 0 0 / 

CALL YRfcAO (NUTZ* V , L V * A V , No T l ) 

return 

end 


WPRT FI.GNATZ 2 



A2-6 


SlKlSHFQOGOO*Fl*GNXY22 

1 SUBROUTINE GNX*22 

2 C 

3 C MAIN PROGRAM TO GENERATE (XU), (JUOF1, ANO (EUL)i AND STORE MATRJ C £S 

9 C ON UTILITY TAPE, 

5 C DEVELOPED by w BENFIELO, C BODLeY, R PHILjPPUS, R WOHLEN. JULY i V 7 3 • 

6 C LAST REVISION &Y R A PHI(,IPPUS, JUNE 197<(. 

7 C 


b C 

9 C INPUT DATA READ IN this PROGRAM. 


10 

C 


10PT 

FORMAT 

{kb) 

1 1 

c 


IF UOPT *EQ, 6HBA5K ) 

CALL READ <XyZ,NJ»3) 


i z 

c 



CALL RE AO l M <J0UF,NJ,6l 


1 3 

c 


• 

CALL READ ( E o L , N J , 3 ) 


i H 

c 


* 

CALL Fi'»EL ISeE subrt for 

Input ) 

lb 

c 


IF < I OP T .EO, 6HXYZEUL) 

call daYgen (See subrt for 

INPUT ) 

16 

c 


keturn 



I 7 

c 





lb 



COMMON / R»»TaPS / NUTEL 

,NUTXYZ,NuTLT,NuTST,NUTMX,NUTnX*NUTBX 

19 



DATA NIT, NOT / 5 , 6 / 



20 

c 





21 


1 COI 

FORMAT ( | 3 A 6 | A 2 ) 



22 

c 





23 



REVKInO nUTEL 



2*4 



RE AO (NiT,lC0l) I9PT 



2b 



IF ( I OP T •EU.6HBAS1C 1 

GO 70 10 


26 



IF (lOPT .Eg# 6HFLAT p) 

GO TO 2m 


27 



IF UOPT .EU.6HCNT bM) 

GO TO 30 


28 



IF ( 1 OP T .E^.bHCYLCOR) 

GO TO 93 


2 V 



IF (lOPT • E^ , 6HC I RC YL ) 

GO TO 50 


30 



IF UOPT . Etv .6HNASTRN ) 

GO TO 60 


31 



IF UOPT .EQ.6HSNAP ) 

GO TO 7 D 


32 



if uopt .ew.^hxyzeulj 

GO TO fi 0 

' 

33 





NERRORa i 

3*4 



GO TO ’99 



3S 

c 





36 


1C 

call phasic 



37 



rewind nUTEL 



38 



REWIND nUTXYZ 



39 



return 



90 

c 







20 



NERRO R ■ 2 

92 



go To 999 



93 


30 



NERROR 03 

9 H 



Go TO 999 

OKitriNAL PAGE IS 


95 


90 

NERR0R»9 

96 

97 


50 

GO TO 999 

OF POOR QUALITY 

N£RROR=S 

98 



GO TO 999 


99 


60 



NE RRQR * 6 

5C 



GO TO 999 



51 


70 



NERR0R»7 

52 



GO TO 999 



53 


80 

CALL XYZEU2 



59 



Rewind nutEl 



5b 



REWIND nUTXYZ 



b6 



RETURN 





b 7 C 

58 9*9 CALL Zl ti 0M8 » NERROR ) 

5v £no 


WP K T F t *MLCHEQ 
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* J *1 5HPGGQQQ*Fl .MECHEL 


1 

z 

3 

4 

5 

6 
7 
tt 
9 

10 
1 1 
1 2 
1 3 

1 4 
IS 
i 6 
17 
1 8 

19 

20 
21 

22 

23 

24 

25 
2b 
21 

2h 

29 

30 

31 

32 

33 

34 
3b 

36 

37 

38 

39 

40 

41 


SUbftQUTiNE MEC*EG 

COMMON / 8 TAPED / NUTKO,NuTUo.NUTOD*NUTU ,NUT|2 # NUU3»NUT14 
DIMENSION V ( 4000) ,LV ( 4C2G ) #T 1150° » 6 1 » T F ( 600, 6)» A < 405,3 ) , 
*0D(4Q5*6) . JV£C<6) f E 1425.3) .Rl6»6),RTU f 6),EM<6.6> 

DIMENSION TT (6.600) # XR(3>*RV(6,6) ,1 VEC(o) , 1TVC(6C0) 

DATA KV »*R »KJ ( Kg ,KX , *BT/ 

* 4030, 6q 0 * 4f?5, 6t 3 » 1500/ 

EQUlvALENCEd, VI , (T(40G1 ) ,LV ) 

C 

c program to calculate effective mass for slosh 
c developed by p a abbott, December" 1974. 
c 

c form rigid body transformation in global 
c 

2001 format 1// » iox»I3hslo$h mass * ,eu * 4 ,sx » bhcp ( x > ■ ,eu*4) 

c 

CALL REaDIX »NR»NC,KUiKX) 

READ(S e |Ol) n3 » NBL A 0 
1 0 I FORMAT ( 1615) 

c 

C N 3 IS THE last node WITH only 3 DOF 
c nblad is the dqf no, at y3, 

00 6 I *1 »N3 
DO 6 0*1.3 

6 JD< I , J)»J*3*i 1-1 ) 

NoN3+ I 

DO 7 I*N.NR 
00 7 0*1 .6 

7 JD ( I ,J)sJ+3«n3+( I ?N ) *6 

CALL WRlHM(u0»NR|6,6MUD0FUL,KJ) 

call readixr ,nf . nc » 1 ,kx i 
call headim< jytc , nrj,mcj , 1 ,kd) 

CALL RBT«HX # xR|JO,JVEC»T.NR # NRT f NCT,Kj t KBT) 

CALL ftRiTE(T # NRT,NCT ( 6HRGL0bL.XBT) 

C 

C TRANSFORM to local system 

call REaoimi jd .nr.ng # kj ,kdj 
Call REaO (E ,NR .NC#KU ,KX j 
do I OO I » 1 » 6 
1 co 1 vec c j > - 1 

L » 1 


42 

43 


M«0 

MI « 0 


44 

4b 

46 

47 

48 

49 
BC 

51 

52 

53 

54 

55 

56 


Do 10 l « 1 » NR 

CALL ZErO { JV£C 1 1 1 6 . 1 ) 

call Z£ko<r ,6 ,6 ,kd 1 

X R OT » 3 

1FU.GT.N3) KR0T«6 
Ksq 

Dq 9 J* I 9 6 

IF (JD d » J ) -LT .0) GO TO 888 
IF (J,LT,4) Xft(j)*&d»J) 

I F ( JO d , J ) .LE.O) GO TQ 9 

IF ( I • <*T 0 N3 A NO , JO (J . JI .LE.NBLAO ) GO TO 9 

K®K+S 

JVEC ( J ) bK 
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57 

5 6 

5? 

60 

M 

62 

6 J 

M 

6 b 
66 

6 7 
68 

6 V 

70 

71 

72 

73 

7 H 
7b 

76 

77 

78 

7? 

80 

61 

82 

63 

8 H 

6b 
66 
87 
66 
6 V 
93 
VI 

92 

93 
VM 
9b 

96 

97 
96 


MI a M i 4 1 

l T V C ( M I 1 ■ I , J) 

9 CONTINUE 

call eulerur,r,kd) 

IF (KR0T.EQ.6) CAUL EUlER ( XR , R (**, H) , KO 1 

call zero { rv ,krgt ,kd ) 

call re*add d« » r,i vec, jyec,rvi*rot ,krot ,krot,k,kd,ko> 
CALL TRANSCRVtWT.KROTiKiKOiKOI 

call mu l ta (rt , t ( l , i ) ,s,krot ,nct ,kd ,kbt » 

Oq 6 IT»I,K 
DQ 8 JT»I ,NCT 
6 TF(M4lT # JT l»KT< 1T,JT 1 

M * M 4 (<, 

86 B L*L+kROT 
1C CONTjWUfc 

call WRJTM Tf ,H»NCT , 6HR0TKAN ,KK 1 

00 889 j*l i3 
$89 OVUC l 0 ) « J 

CALL Zlr 0 ( T » M , NC T , KB T 1 

CALL HEvaDD (l«'TF t ntfCtJVEC»T»H»NCTtH t NCT,KR,KBT) 

CALL WRITE (T »M*NCT,6HC0RTRM t K6Tl 
CALL TRANSIT ,M*NCT (K3T »KO) 

C 

C READ MASS AND MupES AND DO THE KES] 
call yrlad inuti i » v ,lv ,kv ,nuti*o 
CALL TREAD inUT12*V ,LV ,0,NuTH 1 
CALL TMULT1 ( NOT 1 1 *NUT 12 »NUT t 3 » V ,LV ,KV ,NOTiH) 

call tstod inuti3,tf ,nr ,nc ,kr »*d »v ,lv ,kv ,nuti 

CALL MUlT (TTbTF»RiNCT f NR # NC,WO f KRl 

00 lb 1*1, nct 

OQ 15 0“ 1 ,NC 
15 TTC 1 , Ji * H 4 1 ,01 

CALL TKANSm.TF t NCT,NC,K0,KK) 
call w r i r e < t i , ft c t ,nc »6httmphi , k u ) 

Call palEho 
Do 20 1 3 l *NC 

call nult ( ttu , i i ,tf u , n ,em,nct , j ,nct »kd,kr> 

c 20 Call W *1 TE ( Etf ,NCT ,NCT ^AMt IHhMASS,! ) ,KD) 

CpREbS a EM(2 ,3)/£h(2 t 2) 

2c write (6,2001 ) EM( 2, 2) .CPRESS 
RETURN 
END 


feipRT F1#0F1 NEl 
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SiK !SHP0OQO0*f‘| .OF INEL 


SUBROUTINE 0 F l N E L 

MAiN PROGRAM TO READ (XY*), (UOOF), lEUL) ANU CALCULATE (On UPTIQNj 

assembled finite Element mass* stiffness, load transformation, and 
buckling matrices. 

calls foR*a subroutines finel ,yin »ywrite* 

DEVELOPED by B 6ENFIELQ, C BOOLeY, H PHiLjPPUS, R WONLEN* JULY l?;3 

last revision by R a philippus, December 1974. 


10 

1 1 

12 
1 3 


COMMON / K«TApS / nutel,nutxyz,nutlt,noTst,nutmx,nutkx,notba 
common / BTAPEF / NUTMfNUTK ,NUT1 ,NUT2 ,nUT3 

COMMON / DUMMY / V( 12300) 

COMMON / I 0 U h M Y / L V ( 1 2DCC ) 


1 4 

is 

id 
1 7 


COMMON / RESTAP / MRSVT1 

DIMENSION XY2(200P,3>* JuOF ( 2000 , 6 > , EUL(20QC,3> 
EQUIVALENCE (XYZ.V), (EUL.V (6001 > I * (JDOF,LV) 

data kkx, kcx, kru, *cu» kre* kce, kv / 


Id * 2000, 3 t 2000, 6, 2000* 3,12000 / 

19 C READ XYZ»UOOF,EuL FROM nuTXYZ CREATED IN PROGRAM DATGEN. 


20 

Rewind nutxyz 

21 

HEAD 

(NUTXYZ) Nj f NCX.NRJ iNC J ,NRE ,NC£ 

22 




23 

IF (nCX .NE 

• 3 1 GO TO 999 

24 




2 S 

IF (NRU .NE 

. NU .OR. NCU « NE • 6) GO TO 999 

26 




27 

ip (nRE ,n£ 

• NJ .OR. NCE *NE. 3) GO TO 999 

2b 




29 

IF ( 

NJ oGT. 

K R X .Of*. NCX.GT.KuX *0R. 

30 

* 

N R U « 6 T « 

KHJ *0R. NCU.GT.KCU .UR. 

31 

* 

n«e .gt • 

KrE .OR. NCE.GT.KCE) GO TO 999 

32 

head 

(NUTXYZ) ( (UDOF U ,U) , I = i ,NKU) ,U»l ,NCU ) 

33 

READ 

(NuTXTZ) (( XYZ(X,U),I»i, «U).U«l,NCX) 

34 

read 

(NgTXYZ) t( EULU , J) , I«I ,NR£I , J«| ,NC£) 

3S 

call 

WTaPE 

t X Y Z , N U , N C X , 6 h X Y Z ,KRX,NKSVT 1 ) 

36 

CALL 

wTaPE 

( JOOF ,NR J ,NCJ .6HU00F ,KRU # NRSVTl) 

37 

Call 

RTaPE 

( EUL ,NRt ,NCE , 6HtUL~ , KRE , NKS V T l ) 

3b 

call 

y inel 

{XYZ,UD0F»EUL.NUTEL,NU» 

39 

« 


NUTM ,NUTK,NUTLT ,NUTsT,nuTbx, 

40 

• 


V ,LV , K V , K R X ,KRU|KRE , 

4 l 

* 


NUIMX,NuTKX,NUT1»nUT2»NUT3) 

42 

call 

YrtKlTE 

(WUTM,4HMASS , V ,LV iKV ) 

43 

call 

yw T APE 

(NUTM,6HhASS , v , lv ,KV *NHS VT 1) 

44 

CALL 

YffR I TE 

(NUTK,4HbTlF , V , L V , K V ) 

4b 

call 

LTaPE 

( NRS V T l ) 


NEKrUr* I 

NEKRU^»2 


NtHRD K -a 

NEKR0r.4 


46 RETURN 

4 / C 

4b 999 CALL ZZtoOMfa { 6«UF I N£t ,Nt KKOR | 

4 V END 


f«ipR T FI.OPlTIO 
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SiRiSHPQOQQO*F i *OPLT 10 

1 ' subroutine OPLUQ 

2 COMMqN / BTAPEU / NUTfcD »NUTL0 »NUTOD,NUl 1 ,NUTP ,NUT3 ,NwT9 

3 DIMENSION XY](95,3l, Pl 0C<90, 20>. *fM2cO), YP<200)'» aPYP*90,S)» 

4 * DX(95), 0Y<45), DALlRb), OYl|SE»), PT ITLE I 9 ) » 

g * SCAUEMI20)* FKfci<H(2G}» B IGP I 600 , 6 I , I VEC I 95 I , ODOF { 600 9 6 1 • 


6 * V(hOCO), LVi90G0> 

7 DATA *J, KP» *M» KPLOT,KtiIOP ,KBU«» KV/ 

b * 95, 90, 20, 2GQ, 60C, 6« 90 q0/ 

9 data eps/i *e-i&/ 

iq data ni t ,n0t/5»<>/ 

11 l 0 w 1 FORMAT { l 3 A I5»IS) 

\z 2qc i format uo*, f iq*9 »f io.h > 

13 2002 FORMAT (7/lGX, 9HM0DE NO. 13, //» 

14 2o03 FORMAT 1//1QX, 9HM0DE NO. 13, 5* 7HSCALE ■ Fa,3l 

lb 5001 FORMa t I ShMOuE* I 3 » l* 5HFRE0-F? * 4 , IX &HSHFT « l 9 * IX 4HS CL = F 9 • 3 > 

16 CALL JDFRMV (l2H 6 1 flULTMAN , 1 2« BIN 5 1 9Q »*HS 632b , 

17 * I 2HHARDCQPY ) 


lb 

1 9 
2G 

2 1 
22 

23 

24 

25 

26 
27 
2b 
29 

3C 

31 

32 

33 

34 

35 

36 

37 
3b 
3 V 
HO 
HI 


92 

93 
HR 
95 


96 

97 


MB 

99 

SQ 

51 

52 

53 


54 

55 

56 


C 

C READ IVEC. IN SURFACE JO^NT ORDER FROM CENTER OUT* VALUE IS ROW 
C NUMBER ]N aV2 » JOOF iEULER* 

call readim ( jvec# ii,njs, i,kj) 

C READ X Y Z • 

Call READ <bIGP, NJF,I3, KB1G?i3> 

Do S I * l * NO S 
jnt * ivectn 

X Y T < I *n a aiGFMJNT,!* 

5 XYT( 1*21 = BIGP1UNT,21 
C READ EULtR. 

call Read ibisp, njf.ij, kbigp,3) 

DO 6 1 = 1 » N J S 
JNT * I VEC 1 I ) 

6 X Y T 1 i * 3 > a BIG p !uNT,3l 

call ^rite ixyt, njs,3, 3hxyt, kj) 

IMS = 1 
1SHIFT ■ 0 

CALL READ (FREqM, NM»lltXM,l) 

c READ MODES. FIND LARGEST VALUE (ABS) IN EaCH MODE* ISCALEMj, 

CALL YREAO INUTP , V ,LV ,KV ,NDU » 

CALL Y5T0D (NUTP,BlGP,NBIGPfNM,KBlGP,KBIGM,V , L V * K V , N U T | ) 

DO 9 U®1 , N M 

B I CM * 0«0 

DO a 1 = 1 * N B I G P 

VALM = ABS(B1GP( l» J> ) 

IF 1 V ALM ,GT. BIGM) BlGMaVALM 
B C ON T 1 NUE 

9 SCALE** J) * i,/61GM 

CALL WRITE (SCALER, J,NM» 6mSCALEM, 1} 

C HEAD JDOF » 

call readim ijdof, njf» 16, KdiGP,6) 

c compress bigp into ploc. order is u,v in joint order given by ivec. 

DO 12 1 * I | N J s 
JNT a IvECU) 

IU * JDOFI JNT, n 
IV * JD0F(JNT#2) 

Ipu a 2 * 1-1 

JPV ■ 2*1 
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67 DO J 2 3» * j N m 

58 PLOC< IPU, J> a 0*0 

69 IF C I U ,GT* 0) ' PLOCI 1PU,J)«8JGF i iUt J> 

60 ploc<ipv*j) » a®o 

61 IF U v ,GT« Q), PL0CUPV 9 J)*B1GPUV»JI 

62 12 CONTINUE 

63 CALL WRITE <PL0C» 2*NJ5,N«, 9HPL0C. Kp). 

64 C ' 


6b 
6 6 
6 7 
68 

6 V 
7Q 

71 

72 

73 

74 

75 

76 

77 

7 6 

7 V 
80 
81 
82 

S3 

84 

6b 

86 

87 

88 

8 V 
90 

9 I 

92 

93 

94 
9b 
96 
9 7 

98 

99 
*00 
1CI 
*02 
*03 
1 04 
} OS 
1 06 
|07 
108 
109 

no 
* 1 1 
* 12 
i 13 


C DEFINE UNDISTUR6E0 SHAPE® XP,YP IN PLOT X # Y* LEFT i R * GHT • 

NJ o NJS 
NjM* a NU«| 

1 J at N J+ * 

DO IQ *Pol»NjH* 

1 3 ® i o • * 

XPUP* » - X Y r ( 1 0 i 2 ) 

10 YP UP) * XYT(lUsl) 

00 JS IU»l*No 
Ip a NJfil^IU 

XP (Jp ) = X Y T * I J # 2) 

15 YPUP) = XYTU J »* > 

C 

DO 199 IMb1 5 nM 
MODEnO o IMS-1+IM 
SCALE * SCALENUM) 

FREQ « FREQM(IM) 

c CALCULATE DELTA-X# DElTA<?Y, IN JAG X,Y* tEF 7 *DXL., DtL . RIGHT* 
call pagehd 

WRITE <NOT*2Q02) MQDENO 
00 SO I 3“ 1 9 NO 
IXL * 2*U-1 
IYL ® 2*1 3 

ANG = XrT( I J 6 3)/57,2957795l 
CJ b CObiANG) 

SJ « SINIANG) 

DX3AG = PL OC ( lXL v |M)»CJ - PLOC (1 YL » I H ) »$0 
DY3AG PLOCUXL, PLOC (1 YL » I H ) *CU 

WRITE <N0T 9 20Gn dyjag.dxjag 

DXUJ) = XYTuJtl) * 5CALE*DXJAG 
OXL U 3>«*XYT<I3 tl ) - SCAL£*DXJAG 

OYUJ) ® XYT( I U » 2 ) «■ SCALE*DYJA& 

50 OYL( 1 U I * - X Y T ( l 0 * 2 ) ♦ SCAtE*DYJAG 

C PACK DY 6 DX INTO XP»YP« LEF T » R I GhT • 

I P I « NjH*^Nj4l 
IP2 s l P 1 ® 1 +N JM I 
1 3 «a N3-U 
DO 60 I P*IP 1 » I p 2 
1 3 ® I 3« | 

XPUP > ® D Y L 1 i J I 
60 ypjip) a DXL( 13) 

Ip* E> NjMI+NJ+NJM 1 

DO 65 lU» 1 *NJ 
IP a I P * + IU 
XPUP) = D Y ( I J ) 

65 Yp UP I a ox U j) 

C PACK XP » TP FOR PRINTING, 

N J2H * « 2*NJ-j 

DO 60 I B 1 «NJ2M* 
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i lv 


xp yp < 1 1 1 > a xp < i) 

lib 


xpYPd 1 2 1 * yp< n 

i 1 6 


XPYP 1 I * 3 1 « 0*0 

t i 7 


X P Y P < 1 » i> * XPU*Nj2hl) 

& i & 

8C 

Xpyp(I»5) » YPU + N J2M 1 ) 

1 J? 


CALL »«ITE IXPTP, NJ*M1,&, VHXPYP, KP ) 

120 


WRITE IN0T»20C3) MODLf^O , SC ALL 

1 2 1 

C PUT 

IN TAN N TOP , fiO T AT CENTERLINE* 

1 22 


XP( 1 p+1 ) « o* 

123 


Ypllp+l) “16* 

IZH 


Xp { I p ♦ 2 

12b 


T P i I p + 2 l « 0* 

126 


NRP « l P+2 

127 


iFSAflE * I 

128 


IFCORV * i 

12V 


IPLIFT = 1 

130 


ENCODE (HQy&ocl *PTlTtEl MODENQ , F REti , I Sh I F T * 

131 


CALL PLOT 1 | XP * YP *NRP , 1 » -10**2*, 6HY-TANK 

132 

* 

iFSAME.iPCURV, IPLIFT ,KPLGT 1 

133 

19? 

CONT 1 NU£ 

l 3 S 

C 


135 


call pltnd 

136 


return 

137 


end 


^PTjiTUE 


tfPRT Fl.OrMOOE 
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SiKJShPOOOQQ*Fl *0YM0DE 


1 

2 

3 

4 
b 
6 
7 
B 
9 

i 0 
i i 
12 
1 3 
I 9 
lb 

16 
l 7 

16 

19 

20 
21 

22 

23 

29 

2 b 

26 

27 

26 

29 

30 

31 

32 

33 

39 

3b 

36 

37 

36 

39 

9 l 
9 2 
93 
99 
9b 

96 

97 

98 

99 

&0 
b l 

52 

53 

59 

55 

56 


MARCH 1972* 
02APR79, 


SUBROUTINE QYMODE 

MAIN PROGRAM TO TEST ITERATIVE KATLEI 6H-R I TZ METHOD OF D** OQHN A^M* kE 
technique = composite STRUCTURE. 

VERSION * NQNtSttEEPING, 

PROGRAMMING LOGIC ■ SPARSE, 

MAXIMUM size OF mA5S|STIF o 1820* 

MAXIMUM NU » 70 

DEVELOPED BY R L WOHLEN *ND R A PHItlPPUS. 

LAST OJmENSION CHANGE < FOR 1500001 6Y HRuq A 
LAST REVISION BY R A PHILIPPUS* DULY 197m* 


INPUT DATA READ IN THIS PROGRAM, 

N«| FORMAT U0A,I5) 

NU FORMAT CiOX.ib) 

'SHIFT FORMAT ( 1 0 X , El 01 

M A X I T FORMAT ( 10X , 15) 

1 punch FORMAT ( A 6 1 

RETURN 

DEFINITION of INPUT variables. 

Nfl S NUMBER OF MODES WANTED, 

nu « number qf rayle igh-r i tz mooes to use, 

SHIFT « SHIFT VALUE TO USE* 

MAMT - MAXIMUM number of ITERATIONS Tg Qg PERFORMED* 

IPUNCH b PUNCH CARD OPTiQN FOR *V2, FRE$, MODES , AND TMOOES, 

« 6HPUNCH • FOR PUNCH CARD OUTPUT* 
a 6HN0PUNC» NO PUNCH CARD OUTPUT, 

COMMON /BTAPEM/ NU7M,M0TK,NUTTR»NUTZ,NuTF,NUT|»NUT2,NUT3*NUT9 # 
» NUT5 * NUT6 , NUT 7 

COMMON / RESTAP / NRSVTl 
COMMON / RTRaNS / I F TRAN 

DIMENSION V < 1 0920 1 > LVU09201* *2<7Q), w<7Q)» FREQI701* 

« A l 70* 701, S( 70, 701, Mhllo) 

EQUIVALENCE ( V ( 36 9 | ) ,$) , <LVl369l),AI 

DATA NIT, NOT / 5,6 / 
data kv, ka / 

* 10920, 70 / 

data NITERI, N1TER2, 

9 0 , 1 , 

Data IFpRNT/ioOO/ 

lOOl FORMAT UOX, 915) 

1010 FORMAT i 10X, ElO.pl 
l o 2o format 112 A 6 ) 


TOLZ, T0Lft2/ 

1,E-Q6, 1.E-Q9/ 


read 

read 

Read 

read 

read 


INIT, loci ) 

(NIT ,10UI 1 

(NJT.IOIO) 
(NIT, 1001 1 
(NIT ,1020) 


Ntt 

NU 

SHIFT 
MAX I T 
1PUNCH 
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5 7 

58 

59 
6C 

6 l 
62 
63 
6 H 

65 
6 6 
67 
6 6 

6 V 

70 

7 1 

72 

73 
79 
7b 

76 

77 
76 
79 
00 
B1 
02 

03 

89 

05 

06 
07 
06 

89 

90 

u 

92 

93 
99 


CALL YMQDE 2 ( NUTM , NO T K ,NUT 2 , « 2 , W »FREW , N# » V,LV»AiS, KV.KA f 

* N<JTl f NUT2,NUT3,NUT9 > NUT& tf vjUT6 # NUT7» 

* ifprnt,haxit, 

* NO *N I TERJ |NJT£R 2 »SHJFT *T 0 lZ,T 0 LW 2 ) 

CALL WRjT£ < * 2 ,NW , 1 , 2 H« 2 »KA > 

CALL WRITE 16 iNW.IpiHW ,KA) 

Call WRITE (FREQ »NA , I t^HFREy |KA ) 

Call YWRITE < NUT 2 , BHhQDES , V ,LV *ICV ) 

IF (nRSVTI *GT« 01 CALL wTAPE (« 2 tN*» | » 2 H* 2 »KA 9 NRSVT 1 > 

IF (isRSVTi # G T • 0 ) CALL WTAPE ( FREQ # n« , 1 , HHFREQ » K A ,nRS y T 1 ) 

IF InRSvTI ,GT. O) CALL YATAPE tNUT 2 ^HM 0 DES t y |LV ,KV ,N«SVT II 

call ydtos ( fReq ,nut f ,nw , i ,ka,i , v ,lv ,kv ,nuT 7 > 

IF (jFTRAN .EQ* 0 ) GO TO ICO 
REWInO nUTTR 

CALL YlNl (NOTTR i H h # I , 10 ) 

NRT « Mh ( l i 
NCT ■ Hh< 2 > 

CALL. YOJSA (N 0 T 2 . i , 1 t NUTl *NCT ,NW, V ,LV ,KV,NUT 7 ) 

CALL VNULT ( NOTTR ,NUU ,NUT 2 ,V ,L 9 ,KV ,NUT 7 ) 

RE W IND NOT l 

CALL YINI (NOTZ#NRM, l , | ) 

NR TM 8 NRii - NCT f NR | 

call yzerq (nutj *nrtm,n«ij 

Call YASSER (N 0 T 2 , l » ) #NUU • V ,LV iKV»NUTS |NUT 6 ,NUT 7 ) 

NCT P 1 8 NCT ♦ I 

NRX = N^M - NC T 

CALL YDJSA (NOTZ.NCTei ,1 , NU T 3 » NRX » N« » y , L V , K V , NUT 7 I 
NRTP1 8 NRT ♦ I 

CALL YAbSEfl (N 0 T 3 ,NRTP 1 »I *NOTl » V ,LV ,KV f N 0 T 5 ,NUT 6 ,NUT 7 ) 

call y write (nvti« 6 Htiiodes 9 v»lv»kv) 

IF (nRSVTI .ot* o> call ywtape ( NUT 1 , 6 hTM 0 DES 9 V ,LV ,KV ) 
ico if (nRSvti « g t ♦ o) call ltape (nrsvtii 

IF UPOnCH 9 nE • 5 HPUNCH ) RETORN 

call punch (FReq,nw,i *hht re^.kai 

CALL YPUNCH ( NUTZ » 5 HM 0 DES » V , L V ,R V J 

if < iftran »ne* 0 ) call ypunch (nut i , 6 hTmooes , v *lv »kv j 

RETURN 


95 C 

96 END 


®PRT F 1 .PBASIc 



-0 tr 


A2-I6 


5 I NiShP0D0GQ*F 1 • Pd ASK 
1 


10 

u 


SUaRoVTlNt PHASIC 

program to call basic, 
COMMON / RftTApS 


t 

/ r\«iftr j i NUTEL|NUTAY2fNUTLT,NoTST»NUTNX,NUTKX>NOTBX 

DIMENSION X Y Z ( 2000 » 3 ) , JDOF (2000, 6 » I EuL(2000*3) 

Data KRx. KC*. kRj. kc j* krE, K.CE / 


data k«x, 

« 2 CCC 1 
CALL BASIC 

RETURN 

END 


KC*» kRj, kco, 

3 , 2000 , 6 , 2000 , 3 / 

(XYZ, JOOf ,EUL,KRA,KCX,KRJ,KCJ,KRE,KCt;, 
NUTLL ,NUTXYZ ) 


WPRT FloPXTRA 
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* I* IShPOOQOO*F 1 .PXTRA 

1 ' SUBROUTINE PATHa 

overlay program TO R £ V A o o extra mass and STIFFNESS matrices to 
existing mass and stiffness matrices. 

input data read in this program. 

CALL REaDIM (lJVEC.i.NCII 

call tread (mass matrix) 

CALL TREAu (STIF MATRIX) 

RETURN 


2 

3 

*4 

5 

6 
7 

a 

9 

10 

1 1 
l 2 
1 3 
l *4 

15 

16 
1 7 



»*• 


# 

READ 


” s. i v n 

COMMON / BTAPET / NUT M I NUTK I nut l # NUT2 t NUT3 f NUT^ TX 

dimension v<i20Q0)» lv<i2000), uvec(2ooq) 
data kv» kIv / 


I 2030 > 

1 I v; v C F n R 


2C 00 


I 6 
l 9 
20 

22 

23 

2*4 

25 

26 
27 
26 
29 
3 C 

31 

32 


DO 10 1*1, NCI 

IF UUVECCU * NE • 0) GO TO 20 
10 CONTINUE 

DO IS 1*1 • NC I 
lb IjyEC(I) * I 
C READ EXTKA MASS MATRIX. 

20 call tread i not* , v , lv »k v » nu m 

Call YREVAO c 1 • »NUTX , J JVEC , l jVEL ,NUTMt V »L V ,NUT 1 »NuT2 t NUT3 »nD|*) ) 
C READ E X T N A STIFFNESS MATRIX. 

CALL TREAD iNUTX»V,L¥»KV,NUTI) 

CALL TREVAO { | « »NUTX, IJVEC ilJVEL » N U T Ps » V »LV »KV » N U T 1 tNwT2»NUT3»NUT 4 4) 

c call ttirite <nutm ,shmxtha ,v,lv;kvi 

C CALL YKRITE (NUTK .SmkXTRA ,V,LV»FV) 

Return 

End 


Id>P R T T l .REDUCE 
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siRishpgqooo^fi 


1 

2 

3 

4 

s 

6 

7 

et 

9 

10 
1 1 
1 1 

13 

14 
lb 
16 
17 
16 
19 

2 C 
21 


C 

c 

t 

c 

.c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 


c 


• REDUCE 

SUBROUTINE REDUCE 

GVlRlaY PROGRAM TO REDUCE STIFFnESS AND Mass MATRIX. 
OOF TO &E REDUCED MUST BE POSITIONED FIRST IN MATRIX. 
REDUCING TRANSFORMAT I ON IS STORED On NUTTr. 

DEVELOPED BY WA bENFIELO, FEBRUARY 1 974. 


* * *■ * flp a 


input data read in this program. 

NR FORMAT UOX.lb) 

CALL REaDIM UUVEC, J ,NCI 1 

return 


DEFINITION of INPUT VARIABLES* 

nr « number of how-cols in reduced matrix. 

ijvec » i j v ec to Rearrange rows and cols before reducing* 

COMMON /BTAPER / NUTMR.NUTKR.NOTTR.NOTRI ,NUT«2,NUTR3 ? NUTRH,tMUTRs* 
• NUTR6 

COMMON / RTRANS / iftran 


22 

23 

24 

2 5 

26 
27 
2 b 

29 

30 

31 

32 

33 

34 
3b 

36 

37 

3 & 

39 

40 

41 

42 

43 

44 
4$ 

46 

4 7 
4b 
4 9 

50 

51 

52 

53 

54 
5b 


DIMENSION v< 20000 ! I LVUCCOO) 

DIMENSION I V 1 ( 30qC ) , IV2I3QQG1 
C 

DATA nit. not / 5* 6 ? 

Data kv / ioooo / 

DATA *lv / 30&0 / 

C 

1 OO 1 FORMAT 1 6 ( I OX * I 5 I ) 

2 00 1 FORMAT I/// 20X14HOVERLAY REuUCE, 

* // 10X42MNUMBER OF ROWS AND CQLS BEFORE REDUCING » ,!&, 

* / 1QX42HNUMBER OF ROWS AND COLS AFTER REDUCING ■ ,161 

00 to 

1C I V 2 u > fi I 

c 

C READ NUMBER OF HOWS IN REDUCED MATRIX* 

read in i t » looi i nr 

REWIND NUTKR 

call yini inutkr.wrk , j , u 
call pagehd 

WRITE <NOT. 2001 1 NR K > NR 
C 

c READ REARRANGING IV£C. 

call REaDIM UVI »NRl lN CI .1 ,MV) 

DO 12 I*l»NCi 

IF ( I VI < I ) .NE« 0) GO TO 17 
12 CONTINUE 

DO lb 1*1 . N C i 
lb I V 1 < I > * I 
C 

C REARRANGE mass, STIFF, and TRANS <IF AN Y 1 MATRICES. 

17 call yzero inutri ,nrk ,nr< i 
call yeero inutr 2 ,nrk,nrk> 

CALL YREVaD { 1 .0 ,NUTMR , I V 1 ,1 Vl ,NUTRl ,V , L V , K V , NOTH 2 , NUTK3 , NUTK4 » 

* NUTRS) 


56 CALL YREVAD ( 1 .0 ,NUT*R , I V 1 i I V 1 ,NUTR2 ,V ,LV ,KV ,nUTR3 ,NuTR4 ,NUTR& » 
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57 

0 


nutmr } 

56 

IF (IFTKAN • 

EQ • 01 Go TO 2C 

5? 

rewind nuttr 


50 

call 

YUI 

(NUTTR t LV| 1 i 2 ) 

5i 

call 

YZfcftO 

< NUTR3 j L V A 1 ) *L V ( 2 ) ) 

52 

call 

YREVAO 

< l eO,NUTTR.IV2ilVl *NUTR3 ( V f LV l KV,NUTRH*NuTR&iNUTMR* 

53 

* 


NUTKR ) 

6H 

20 CALL 

YSkE02 

(NUTR2*NUTKR»NUT1R,NR,1 *V ,lV ,KV ,NUTRM > NU lR5 f NUTR6, 

55 

* 


NUTMR ) 

56 

CALL 

Y W R I T E 

(NUTTR,2HTR ,v*lv,kv> 

57 

CALL 

YWkITE 

(NUTkR,2hwR iY»LV,KV) 

56 

IF ( i f Tran • 

EG* 0) GO TO 3C 

59 

CALL 

YHULTB 

(NUTR3»NU| : TR|V » L V jKY iNUTRM^NUTRS) 

70 

CALL 

YWkITE 

i nuttr »6HMULTTR *v ,lv ,KV) 

7 l 

30 call 

YBt AB 

< NUTR1 #NUTTR , nutmr >y »LV »WV ^UTRH^UTR&I 

72 

CALL 

Y Ar I TE 

(NUTMR»2HHR p V*LV,KV) 

73 

IFTRaN « I 


7H 

RETURN 


7b 

E.ND 




WPKT F 1 * S a V L 
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51R15hP00QG0*FI*SAVE 


1 

2 
3 

*1 

b 

6 

7 

b 

V 

1C 
1 1 
1 2 
i 3 
IH 
lb 
1 6 
l 7 

U 

1 9 
20 
21 
22 
23 

2N 

2b 

26 

2 7 

26 

29 


SUBROUTINE SAVE 

OVERLAY PROGRAM TO SAVE. MASS, STIFFNESS, aND REuUCIN 0 Tr ANSFOKM AT 1 0* 

MATRICES ON FORMA LIBRARY TAPE. 

DEVELOPED 8 Y WA SkNFlELD* FEBRUARY 197*1* 

INPUT PATA READ JN THIS PROGRAM* 

NaMEN»NaMEK,NA m ET 

return 

definition qf input variables. 


» « « * 9 * * 


FORMAT (3A6) 


N A M E H 


NAM Ex 


NAMET 


name of 
6H 

name of 

4>H 

name of 


MASS MATRIX TO uSt ON FOrMA TAPE* 

» MASS IS NUT WRITTEN UN FORMA TAPE* 

STIF MATRIX to use pN FORMA tape* 

» STIF IS NOT WRITTEN UN FORMA TAPE* 
transformation matrix to use on forma tape* 

, TRANSFORMATION is NO 1 WRITTEN on form tape* 


common / b Tapes / nut^5»nutks» nutts 

common / restap / nksvti 

COMMON / KTRanS / IFTRAN 

DIMENSION VUSOOOF, LVilbCOO) 


DATA *V / 

data nit » 


isuoo / 

NOT 7 b , 


6 / 


ICtH FORMAT ( 1 2 A 6 ) 


3C 
3 1 

32 

33 
3*1 
3b 

36 

37 
3b 


READ 1 N I T , 1 C g 1 ) NAMEM,NAIlEK,NAMtT 

IF UjAMfcn • Nfc. * ) CALL YWTaPE 1 N U T H S , N A* E« » V »L V . R V , NRS V » 1 1 

IF INAMtK *N£ • 6H 1 CALL Y«TAPE 1 NUTK S , N A*EK , V , L V . K V , NRS V T j, 1 

IF l iFTkAN • nE • 0 • A NO * NAME! *NE* 6H ) CALL YWTAPt |NUTTS t 

NAMET ,V ,LV ,KV ,NRSV T 1 ) 

CALL LTaPE (NRSVT 1 ) 

RETURN 

end 


WPRT FI * SUb 5 Tr 
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$>1*1 SnPQOOGQ * F l .5U&STR 


2 

3 

*4 

5 

6 
7 

a 

v 

10 
1 1 
12 
1 3 
I 9 
lb 
16 
1 7 

1 H 

1 9 

2 C 
21 
22 
23 
29 
2b 
26 
2 ; 
26 


SUBROUTINE SUBSTR 

OVERLAY PROGRAM to COMBINE SUBSTRUCTURE MASS AND STIFFNESS MATRICES 
TOGETHER* 

DEVELOPED 0'r WA aENFICLD* FEBRUARY 1 979 . 


!«•«««« 


FORMAT 

FORMAT 

FORMAT 


( )0X, 15) 
UOX.IS) 
( A6 1 


FORMAT { 1 0 X t 1 5 ) 


*»***»»»♦ 

INPUT DATA READ JN THIS PROGRAM. 

NSU8S 
NR 

NTRAN 

00 30 1 * NSU B 5 

NSRC 

CALL YREAD {MASS MAT R I X 1 
30 CALL YREAD (STIF MATRIX) 

IF (NTRAN .EU. 6HTRANS ) CALL YKEaD (TRANS MATRIX) 
return 

definition of Input variables* 

NSUBS - NUMBER O f SUBSTRUCTURES to be kE ad in* 

NR = NUMBER oF ROW-COLS IN TOTAL MASS-STIF MATRICES, 

NTRAN » OPTION TO READ IN TRANSFORMATION MATRIX* 

s 6 h T R A N S » TRANSFORMATION MATRIX JS READ IN* 

= shnotrani transformation matrix js not read in. 
nswc = start how-col to assemble substructure into. 

common / BTAPEP / NUTMP,NUTKP,NUTTP*NU T iP > NUT2P*NUT3P,NUTHP,NUTsP 
common / RTRanS / iftran 


2 V 

30 

31 

32 

33 
3*4 

35 

36 

37 
3b 
39 
90 
9 I 


C 

dimension vuooao)» lvugoco 

DIMENSION I V 1 ( 30C0 ) 

c 

data nit. not / 5» 6 / 
data XV / 10000 / 

DATA K1 V / 3000 / 

c 

ICO! f ORMaT (6 ( lox , IS) 1 

1002 FORMAT i I 2 A 6 ) 

c 

C READ NUMBER of SUBSTRUCTURES TO COMBINE TOGETHER. 
READ (NlTiiCOD NSU&S 


92 

93 
9*4 

95 

96 

97 
9 b 
99 
50 
S l 
5 2 
53 
59 

55 

56 


C 

C READ NUMBER Of HO w "COLS IN FINAL MATRIX, A ND IF TRANSFOKMAT ION WJ(*L 

c be input* 

READ (NIT , IQQI 1 NR 
READ (NIT , I QQ2 1 NTRAN 
call yzero (NUTMP.NR.NR) 
call YZERO (NUTXP ,nr ,nr ) 

DO 30 K = I , N S U B S 

c 

C READ STARTING ROW-COL NUMBER TO ASSEMBLE SUBSTRUCTURE InTO* 

READ (NIT, lOQl ) N5RC 

IF (nSUbS *Ey. I • AND , NSRC *EU* I) GO TO 90 
DO 10 1*1 »NR 
io i v i ( i ) * nSrc ♦ i - 
c 


1 
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5 7 C READ SUBSTRUCTURE HASS MATRIX* 

Stf CALL YREAD ( NUT | P , V ,UV |KV »NUT2H ) 

SV CALL YRLVAL (J *Q,NUT IP , 1 l i 1 V 1 .NUTMP , V t LV ,*V , nUT 2P »NUT3p ,fcUTRP * 

6c • NUT&P> 

61 C 

62 C READ SUBSTRUCTURE STIFF MATRIX* 

63 CALL YRtAO ( NUT 1 P , V , L V * K V , NO T2? I 

65 CALL YRLVAU ( J . 0 , NUT \ P , l V l » 1 V 1 » NUTKP , V t L V , K V , NUT 2P , NuT 3r # NUTHP » 

6 b •» NUT5P) 

66 3Q CONTINUE 

67 CCCCC Call YWrITE {NUTKP ,6HTOTSTF , V *lv ( KV I 

6b CCCCC Call YWrITE (NUTMP ,6hT 0TMAS » V |LV t KV > 

6 V 00 TO SC 

73 CALL YRtAU i NUTMP , V , LV , K V , NUT2P ) 

71 Call YREAD i nUTkP , V # LV ,KV ,NUT2P I 

72 C 

73 C READ TRANSFORMATION MATRIX (IF ANY). 

7 H so IF ( N T R AN • ME * 6 HTRAM 5 I RETURN 

7 S Call YREAD (NUTTPtV,LV*KV.NuT2P) 

76 IFTRAN « i 

77 RETURN 

7 a end 


QPRT M*XYZEU2 
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S 1 H i ShPOOOOQ*F i • 
1 


i 0 
1 1 

12 

n 


XY2EU2 

SUBROUTINE XYZEU2 


SUBROUTINE X Y Z U A(i <NUTEL»WUTAY2 ( 

XYZ< 1200,3 ) ,EUL( 1200 ,3 > ,UDOF< 1200 ,6 ) 

ICONF t 1500 B a I , ICONS (2 00,** ) i l CONB l 2C0 » H ) 

XY23(2aC,3) »EULB(2Q0»3 I 

Kb (20) |M(20) ,X0(2O> #NP I NTS (20) ,N£L£MI (20) **01120) 
AN6( JO) * AN6 i 1 I ) , JNTG(2) 

COMMON / RttTAPS / NUTEL *NUTXYZ 


DIMENSION 

dimension 

DIMENSION 
DIMENSION 
0 1 MEnS I ON 


.NUTLT.fciiTWT 


IH 

15 

16 
17 
) to 

1 v 


DATA PEnTA »HEXA/5HP£nTA tHHHEAA/ 

DATA NCa *NC J/3 ,6/ 

DATA AN6/6HR£TURN ,6HFUU ID |6rt<aKAVTY * 6H , ^HTKNOE ,6H«UA0 . 

* 6Hm 1 t 6HM* 1 6 hK 1 »6H*2 ,feHCN5T /, 

* AN5/5HR0 iSHBeKM .SHGvX ,5HGVY ,5HGVZ »5HE »bhNU 

* 5HTMAS ' ,5HTMEM ,5HTBEN /, 


2D 
21 
22 
23 
2 M 
25 


* I N T o/C # 0/ 

DATA ZZtR0/O.O/ 

data eps/i sct-Ob/ 

C 

C 

loci FORMAT (1615) 


26 1002 FORMAT (QF1Q«C> 

27 1 C'O 3 FORMAT (BAlC) 

2b 3cd format 4 a 6 » 

2 9 3002 FORMAT (5<A6 # HXn 

3c 3003 Format < 3 ( a& # £ io.3 1 ) 

31 3oor format (<415, 3£io, 3i 

32 3C0& Format (515,3ElO*3) 

33 * 3 Q3 6 FORMAT (915) 


I 


3h 

5 oC 1 

F ORMaT 

( ICi / ) 





35 


* 1 bX » 

22HRAD1US 

UF THE SPHERE 

K | 

F 1C • 3 * // 

36 


0 1 5* » 

22H 

FLUID 

MASS DENSITY 

S | 

F I 0 * 7 t // 

37 


* lb* » 

22h 


BULK MODULUS 

V 1 

F i H • 3 | / / 

3b 


• 15*» 

2 2H 


GV* 

31 f 

F 1 0 . 3 » / / 

39 


• 1 5* » 

22 h 


GVY 

9 9 

FIG. 3.// 

HO 


* 15* » 

22 h 


GVZ 

m | 

F 10*3,// 

H 1 


* IS* » 

22HtoEADDER 

MASS DENSITY 

a t 

F 13.7 *// 

H2 


« 15X » 

2 2 H 


E 

n 1 

F1M.3,// 

H 3 


* 1 5 x , 

2 2 h 


NU 

« l 

FJ0.7 ,// 

H H 


♦ lb*» 

2 2 h 

bladder THICKNESS 

a t 

F 1 0 » 7 ) 

H5 

5002 

FORMAT 

(1CM/) 





H 6 


* l b* » 

2 2rt 


ncik 

a t 

lb,// 

H7 


V 15** 

22H 


nkba? 

s | 

lb) 

Hb 

U a 

C 

C 







50 

L 1 

c 

c 


T ^ 





52 

c 

J N r u 

1 v » 1 

- 




53 

t 

READ <5 

» 1 CC2 > 

R 



(BFlO.G) 

5 H 

c 

Read (5 

* 1002 ) 

KHOF , 

BLKM 1 CVX » CV Y * (a V2 

( to F 1 0 • 0 > 

55 

c 

head (5 

* l C C 2 ) 

RHOaLD,E&LD»ANUBLD 

tTBLD (8FJ0.0) 

56 

c 

read (5 

» 1G01 ) 

NCI R t 

NRBAY 


l l 6 I b ) 


« 
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5? 

C 

call read iai ,nrainca,i ,kb) 

(S OF THEM) 


58 

C 

CALL HEAD (Rib ,NRR»NCR , I »K6) 

(1 X N R B A Y ) 


5 V 

c 

CALL KEA0IK (NP1NTS, NRN.nCn, l »K6 1 

(1 X (NHBAY+iU 

60 

c 

call REA01M (NEOEMI ,NR£ ,NCE , 1 1 

{ i X NR8AY ) 


61 

c 

real (s b icC 2> ixom *j-i,npi> 

(OFiO.O) 


6 2 

t 

RE A0 15,1002 ) 1X0X U ) , l»l ,NP2) 

(bF 1C.0) 


63 

c 

00 ZgG l 2»2 ,nr8AY 



6 4 

c 

REAO i 5 » 1 002 ) UOl(K) B H"I »NP1 ) 

(6F J 0.0 ) 


6b 

c 2C0 

continue 



66 

C 

read (5,3002) ELNIO 1 ,ELMI02 

(5IA6.4X) ) 


6 i 

0 

HEAD ( 5 , 3 0 D 2 ) ELM1U1 .ELH102 

1 5 ( A 6 » 4 X ) ) 


66 

C 

00 5q 0 1 l«3,NRBAYl 



6V 

c 

READ (5,3QQ2) ELHlOi ,ELMI02 

IS I A6 t HX) ) 


70 

C BOO 

CONTINUE 



7 1 

c 




11 

c 



i 

73 

c 




74 

c 

E X P L A N a T IONS ---------- 

•• 


7b 

c 




76 

o 

AisCOEFPiC JENTS QF THE FREt- SURFACE POLYNOMIAL, 


77 

c- 

R»R AO 1 US OF THE SPHERE • 



7b 

c- 

ncIR-nq*qf circumferential p 3 5 • 

AROUND the model. 


7 V 

c- 

NRbAYbNOo radial bays in the model. 


3 0 

c- 

rshraojal say radii 



ei 

c- 

THETA1«T0TAL MODEL ANGLE (=180,0) 



62 

C“ 

NPQ I N TSsNO « Of POINTS on EACH VERTICAL LINE IN THE PLANE, 

IVECTQK J 

63 

c- 

XQ-XQ-S pOR FOR THE CENTERLINE* 


84 

c- 

XOl raX-VALUtS IN each vertical 

LINE IN THE PLANE (NO. 

*NP 0 I NT £ ) 

6b 

c- 

(NOTEtPOJNTS ON SPHERE AND 

bladder are replaced 

BY ACTU a u 

66 

c- 

CALCULATED VALVES). STARTING FROM IN EXCEPTING 

cent,li n e. 

67 

c- 

ELMJDIbFIRST ELEMENTIFROM BOTTOM) 

FOR THE ANNULAR SPACE 

ID, 

6b 

c- 

ELHlD2=LAST ELEMENTiFROM BOTTOM) 

FOR THE ANNULAR SPACE 

ID, 

69 

c 

PEnTAssPENTahEDRON 



Vo 

c 

H£XA=HEXAHE°PON 



91 

c- 

NELEMIaNO.OF ELEMENTS IN EACH DAY 

(PLaNEH MODEL) (VECTOR) 

92 

c 




93 

c 




94 

c- 

BLADDER is SUPPOSED to BE ATTACHED 

AT DIAMETRAL PLANE 


9b 

c 




96 

c 




97 

c- 

KhOF^MASS DENSITY of THE FLUID 



9a 

c- 

■ BlKMbBULk MODULUS FOR THE FLUID 



99 

c- 

GVXsGRAviTT acceleration in X- 

DIRECTION 


1 00 

C “ 

OVYaGRAVlTV acceleration in y- 

DIRECTION 


101 

t- 

GV2oGRA V I T Y ACCELERATION IN 2- 

direction 


|02 

c « 

RHOBED*MASS DENSITY QF THE dL ADDER 


103 

C - 

EBLDsYOUNg-S modulus for the BLADDER, 


104 

c- 

ANUeLDaPuJ SSO^-S RATIO FOR THE BLADDER 


105 

C™ 

TBEDaTHlCKNESS of the blaudek 



106 

c 




1 07 

C-nOT t 

-THE TOP ELEMENTS IN THE FIRST ANNULAR SPACE HAVE TO BE 

TETRAS* 

106 
1 0 9 

c 

c 




1 10 

c 

"* - - 



1 1 1 


REWIND nutel 



1 1 2 


rewind nUTXYZ 



113. 

c 
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1 1 9 
115 

! U 
117 

1 1 6 
119 

I 20 
121 
122 
123 
1 29 
125 
1 26 
127 

1 28 
1 29 
13C 

131 

132 

133 
139 
13b 

136 

137 
136 
1 39 

190 

191 

192 
1 93 
199 

195 

196 

197 
196 
1 99 

150 

151 

152 

153 
159 

155 

156 

157 
156 
159 
1 6$ 
161 
162 
163 
169 

165 

166 
167 
166 

169 

170 


PI « ATAN2(Ca0»°i«0) 

RADD£G * ATAN 2 U*P» 1 * 0 ) / 95,0 
C 

C- INITIALIZE SQMfc VARIABLES REQUIRED 
C 

NOP » 0 
NO = 0 
NEL * C 

NELSs 0 
NElBr 0 
NOPBs 0 
C 

call zero <xyz 9 Kx 9 3 ,kx) 

call 2 EK 0 (XYZa»KS,J,KS) 

CALL 2 ER 0 (EuL »AX , 3 ,K£ ) 

CALL ZERO (EULB »KS, 3 ,KS) 

00 20 1 * 1 t )500 
DO 2 Q 3»| , 6 

20 1 C 0 NE < I # 3 ) ■ 0 
DO 30 l?i>KS 
00 30 3 * 1,6 
ICONS 1 1 P 3 1 s 0 
30 1 CQN 6 l I s 3 ) * 0 
C 

c- neao in oat a 
c 

READ 15,10021 R 

READ 15 , 1002 ) RHOF , BLKH » G V X , G V Y » V Z 
READ 15 , 1002 ) RHO 8 LO,EBLO,ANU 6 L 0 9 T 8 LO 

c 

READ ( 5 , 1001 ) NC 1 R,NRBAY 

c 

WRITE (6I&00I1 RjRHOP ,BLXN,{iVX,&VY f ^VZ,RHOBLO,£bLD,ANUbLD,TbLD 
WRITE ( 6 1 50 Q 2 ) NCIR»nRBAY 
C 

CALL R£ AO (AI ,NRA ,NCA , I ,KB ) 

call read < r q ,nrr ,ncr , i ,kb ) 
call REaojm inpinTs, nrn ,ncm , i ,xb ) 

CALL REaOIM (NELEMI ,NRE,NCE ,1 ,X 8 ) 

c 

N C I R 1 s NClR - 1 
ANCIR 3 NC 1 R 

ANCIR 1 ■ NC I R 1 
THETA 7 ■ PI 

THETaD s ThETAT / ANC INI 
NOV 2 ■ NC 1 R 1 
NO » NC 1 R 1 / 2 
NQl a Ny 0 2 
NQ2 ® Ny * S 
C 

c TOTAL NO* OF FLU 1 P POINTS* 

NTPF = 0 
00 15 K®2,NCn 
15 NTPF. a NTPF ♦ NPjNTSlp;) 

NTPF « (NTPFoNClR) ♦ NPINTS(l) 
c- points in The first annular space 
NP I » NP SNTSU > 
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i 7 I 

172 

173 
\7H 
l 7b 

176 

177 
176 
179 
l&O 
1»1 
IB 2 
l B 3 
l Bq 
1 B b 
1 B 6 
167 
1 86 
169 

190 

I V 1 
l 92 
193 

1 9 «i 

195 

196 

197 
196 

199 

200 
201 
202 
203 
20H 
20b 
206 
207 
206 
209 

2 10 
211 
212 

213 
21 M 

21b 
216 
2 i 7 
2 l 6 
219 
2 2C 
221 
2 22 
223 
2 2h 

22b 

226 

227 


C 


C 


C 


C 


C- 


c- 


c 


Np 2 = NP1NTS121 

REAO 15,1002) 1X0UI ,l a l,NPl) 
HEAD 15,1002) lXQltI)#I*l,NP21 


YI « 0*0 

Xbl * R - SQHT 1R**2-Y | **2 > 

1 m 0 * 0 

A1 « A I ( 1 ) ♦ All2>*Yl ♦ A I l 3 ) * Y l * *2 * + Ati5)*Yi**H 

ATI * R + SqkT 1K**2-Yi**2l - *1 

Y2 * RO l 1 ) 

XB2 s R - Sqrt IR*»2-Y2**2 > 

W 2 * 0 • 0 

W2 * AI(1) ♦ All2)*Y2 + Ait3)«Y2**2 ♦ Al<9)*T2**3 ♦ AII&MY2**? 
XT2 x R + SQRT 1R**2“Y2**2) - «2 

XQ U ) * X B 1 
AO tN p l) 3 X T 1 
X01 u i = X 6 2 
X 0 1 1 N P 2 ) e X T 2 


NDPB 3 1 

FIRST POINT op THE MODEL 
NO ■ NO ♦ 1 
JUOF { 1 9 1 ) * NO 

EUL (1.3) * (Pl/2eC) / RADDEG 
ALL OTHERS IN THE FIRST BAY* 

NOP a 1 

oo 60 I®1,NP2 

THETA 3 0,0 

OQ 60 0*1 ,NC1R 

14 DP » NOP ♦ 1 > 

xyzinop.i ) « xonn 

XY2(N0P,2J * RB ( 1 ) • COS ( THE TA 1 

XYZ1nOP» 3) a R 8 1 I ) * S1N1THETA) 

IF (1»NE#NP I NTS ( 2 ) ) 00 TO 6b 

EOL1N0P,!) « aTAn21XY21NOP,3) ,XV21NDP,2) 1 / RaOOEG 

EUL 1 nOP | 2 ) * 0*0 

IF (j *bT. n GO TO b2 

TERM * IR**2 -XYZlN0P,2)« # 2) 

IF ( aBS ( TERM) «LT* EPS) DEG * -90,0 
IF IaBS(TERM) *LT* EPS) GO TO B2 

TERM1 = 1 I *0/2*0) * (R**2 -XYZ<NDP.2)«*2) ** i-0«5) • 

* 1-2.0*XYZ1NDP,2> ) 

* - <A1<2) ♦ 2 o Q # A i (3)«XY2(N0P,2) ♦ 3 • 0* A ll H ) * A T 1 1 NOP , 2 ) * # 2 

* ♦ <t«0*Al(S ) *XYZ 1NDP ,2) ♦ •3 ) 

DEG a ATANITEftNl) 2 KftDDEG 

62 CONTINUE 

EUL(nOP,3) * “90*0 - OEG 


NOPb * nDPB + l 
X Y ZB 1 NDpB ,1) a XYZlNOp,!) 
XY2B(NDp6 e 2) a X Y Z ( NOP , 2 ) 
XY2B(NOpfi,3) » XYZINOP, 3) 
tULB ( NDpB , 1 ) » EUL 1 NOP , 1 ) 

EUtB < NOpB t 2 } 8 EUL 1 Nup , 2 ) 
tULB(NDpB,3) « EUL 1 NOP , 3 ) 
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22 8 

229 

230 

23 i 

232 

233 
23** 

235 

236 

237 
23b 
239 
2 SQ 
2 H i 

2S2 

2S 3 
2SS 
2 SB 


6b CONTINUE 

IF ( l *N£. 1 > 60 TO 66 

fcUL<NOp # l) * aTAn2(XY2(N0P»3| |XYZ(N0Pi2I I / RA&DEG 
£UL(N0P S 3) ® -ASINU2*0*(XYZ(NDP, / 

1 ( S«Q* UYZiNDP » * N.U*XYZ (NDP # 3)a«2 + 

2 S«0eXYZ(N0P»2)»*2l**O*5» / R AQDEG 

IF (3 *GT, NG) 2 ) GO TO 79 
IF ij.Ej.l « OR * U.Efl.NQ2> GO TO 76 
NO * NO 4- 1 

JOOF (NOp t \ ) = NO 

&0 ■ NO ♦ I 
JOOF ( NOp # 3 ) = ND 

GO TO ?9 
76 CONTINUE 

IF t J • £ Q • l ) KD ■ -I 

IF (j «£Q« N02) KD « 3 
No » NO ♦ 1 

JOOF(NDp,K0) = ND 


2 S 6 
2 S 7 
2 S b 

259 

260 

251 

252 

253 

255 
25b 

256 

257 
2 b 8 

259 

260 
261 


Go TQ 79 

66 CONUNU'fc. 

IF <J * G T • NJ2* GO TO 79 
IF (j.Ej.i .QR« J.EQ.NG2I GO TO 75 
NO * N 0 ♦ ( 

JOQFlNDp*n = NO 
ND a NO ♦ l 
JOOF 1 NOp # 2 ) a NO 
NO » NO <*■ 1 
JOOF < NOp * 3 I 3 NO 
GO TO 79 
75 CONTINUE 

IF (J •N£. 1) GO TO 77 
NO * ND M 
JOOF ( NOp t J > s NO 
NO ». NO + 1 


2 62 
263 


JOOF < NOp , 2 ) a NO 
GO TO 79 


269 

265 

266 
267 


77 CONTINUE 

IF * I *NE. NP2> KO « 2 
IF <1 *EW« NP 2 1 K 0 * 3 
NQ a NO ♦ | 


26b JOOF(NOP,XO} = NO 

269 79 CONTINUE 


270 THETA * THETa * THETaO 

271 80 CONTINUE 


272 

273 


IF (1 .EG. 1 .JR. I »EJ. NP2) GO TO }0b2 
NOP = NOP - NC I H 


279 

27S 


N l 3 NOp 

N2 a NOP ♦ NO Y 2 + 2 


276 

277 
27b 

279 

280 
2el 
262 
283 
269 


C 


c 


N 1 » N 1 * I 

N2 = N 2 - l 

JOOF 1*2,1) = -JOOF ( N l , 1 ) 
J 0 0 F ( N 2 » 2 ) a JOOF ( N 1 § 2 ) 

00 81 J a I*nQI 
[4 J s N t + 1 

N2 « N2 - 1 
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2bb 

266 

267 

26 S 

269 

290 

291 

292 

293 

294 
29s 

296 

297 

295 

299 

300 

301 

302 

303 

304 
300 

306 

307 
306 

309 

310 

31 1 

312 

313 
3 1 4 

31b 
3 1 6 

3 1 7 
3 1 B 

3 1 9 

323 

32 1 
322 
3 23 

324 

325 

326 

327 
326 

329 

330 

331 

332 

333 

334 

335 

336 

3 37 
33B 

339 

340 

341 


C 

JOGF <N2, 1 > * -JD0F1N1 , 1 1 
J00FCN2 P 2J » JOQF ( N 1 | 2 ) 

JDQF ( N2 fl 3 1 ■ -OOOFIN1,3) 

61 CONTINUE 
GO TO 1084 
1 Q82 CONTINUE 
C 

NQP a NOP <■ NClN 
N 1 » NOP 

N2 ■ NOP ♦ NOV* ♦ 2 
C 

IJ1 e HI + 1 

N2 » N2 - 1 

JDQF t N2 , U » »UDGF (Nl , 1) 

JDOF (,N2 f 2 ) * »JDQF(N1»21 

c 

D 0 1 0 8 3 J a l » N Q 1 
Nl E Nl H 
N2 a N2 - I 
C 

JO OF ( N2 & 1 ) ns «JDQF INI »l 1 
J00F(N2,2| « «UU0F<N1 # 21 
JDOF ( N2 , 3 1 * UOOF ( N 1 , 3 ) 

1 C® 3 CONTINUE 
C 

IQ&H continue 

c 

NOP a NOP ♦ NCIN 
C “ CENTERLINE POINTS 
NDP = NoP ♦ I 
xyzindp, n » xou + n 

NO « NO + l 
JDQF 1 NDp t 2 ) a NO 

60 continue 

C-- FOR THE TOP POINT (TOUCHING BLADDER ) 

JO OF ( NO P g I 1 a JOOF ( NOP t 2 ) 

JDOF(N0p»2l • 0 
XYZB< 1 9 I ) ® XYZ INDP * 1 1 

EULB ( 1 » 3 1 a (-Pl/2,01 / HADOEG 
C 

EuL ( NDP 9 3 I 8 EULB ( 1 9 34 

c- all other node points starting frgh second bay outwards. 

C ' " 

DO 200 II«2 # NRBAY 
NPI s NPJN15( | !♦! > 

READ i S « 1 002 ) UOI (K) ,|Cal ,Npi 1 
t a fiB t I I ) 

xa a R - SyRT i R**2-Y*92 1 

W ffl o • c 

* * A^D ♦ A I C 2 1 * Y ♦ AI(31*r*®2 * Ai(4l*Y*«3 ♦ AJ(&>*Y**4 
XT a R * SORT lR**2-Y»92 ) - <Y 
XOI ( i > a XB 
XO l ( nP l 1 a XT 
C 

DO 120 I » 1 » N P I 
THETA * o oO 
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342 

343 

344 
34b 

346 

347 
346 

349 

350 

35 1 

352 
3b 3 
354 
3&5 

356 

357 

35a 
3 &9 

360 

361 

362 

363 

364 

365 

366 

367 
366 

369 

370 

37 l 

372 

373 

374 

375 

376 
37 7 
376 

379 

380 
361 

382 

383 

384 

385 
366 
387 

386 
389 
3 V D 

391 

392 

393 

394 
396 
396 
39 7 
396 


DO 140 0= | »NC1* 

NQP 7t NOP ♦ | 

XY2 INDP iD » XOI ( I) 

xyz<nDp, 2) * (?8im*cos( theta) 

*YZ(NDP,3) B r 6 < I I ) *5 id t THET A ) 

IF U *Nt.NPlNTS( ll + \) ) GO TO US 

£Ul(nDP,|-) a ATAnZ Urz 1NDP ,3) »XYZ(N0P»2J| / RA0D£G 

£UL(NDP,21 « 0*0 

IF ( J *GT • 1 ) 60 TO 122 

TERM * <R**2 - XYZ { NOP , 2 1 **2 ) 

IF <a8S(T£*M) *LT « EP^ DEG = -90,0 
IF UBS ( TERM) » LT e EPS) GO TO U2 ' 

TERmi a 11.0/2*0) * 1R**2 -XTZ INDP *2 ) **2 ) *• j-0.5) » 

* U2,0*AYZ (NOP ,23 ) 

* - ( A 1 ( 2 ) ♦ 2#0*AI (3)*XYZ(NOP,2) 4 3 • q* A 1 (4)*aYZ(NDP,2)*#2 

* * 4 « 0* A 1 (5 ) *XYZ<N0P ,2 ) **3 J 
OEG a ATAN<TERH1) / kaodeg 

1 2 2 CQNT1 NOE 

EUL(nDP ? 3) « -90*0 - DEG 


NDPB 3 NDPB ♦ 1 

XYZ£U NDpB , 1 ) = XYZINUP,!) 

XYZB(NOp6,2) a XYZ(NDP,2) 

XYZB(NDpa # 3) a XTZ(N0P,31 
EULB(NUpa t l) a EOL ( NDP , 1 ) 

EUUB (NOpS ,2 ) 3 E U L ( N D P , 2 ) 

EUUB ( NOp B , 3 ) e EUL(N0p f 3) 

1 2 5 CONTINUE 

IF (UNe*!) 00 TO 126 

EUU<NDP,1 ) « aTANZ(XYZ(N0Pi3) »XYZ(NDP,2) 1 / raooeg 

EuL(nOP,3) ■ “ASIN( ( 24C* (XYZ(N0P, 1 )-K) } / 

1 ( 4 0Q* ( *Y2 ( NOP • U -R I **2 + 4 • 0*XYZ ( NOP » 3 ) **2 + 

2 4.0*XYZ (N0P »2)**2)6*0.5) / RADDEG 
IF ( J *GT e Nw 2 > GO TO 139 

IF ( J • E u . 1 * OR • J.EG.N02) GO TO 136 

Np c NO * J 

JDOF(NOp,l> « NO 
NO s NO + { 

ODOF ( NOp , 3 ) a NO 
GO TO 139 
136 CONTINUE 

IF(joEQ»1)KD«1 
IF (0 • £ 0 * N ^ 2 > KD a 3 
NO a ND ♦ l 


OOOF ( NOp » KO ) a NO 
GO TO 139 
U6 CONTINUE 

IF (u »GT» NQ2 ) GO TO 139 

IF (J.E <*•! .OR* J.E0.NW2) Go TO 135 
NQ <k NO * l 

0 DOF ( NOp g 1 ) a NO 
NO a NO ♦ 1 

J00F(NDp,2) * NO 
NQ « ND + | 

J0OF ( NOp , 3 ) a NO 
GO TO 1 39 

i35 continue 
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399 


IF <J sNE, 1) GO TO 137 

90G 


NO o NO + i 

9Q 1 


JOOF ( NDp » | ) = NO 

90 2 


NO » NO l 

903 


JOOF ( NOp e 2 1 « NO 

909 


GO TO 1 3 9 

90S 

137 

continue 

906 


IF ( t «NE # NP 1 > KO « 2 

907 


IF < I • £0* NPI > KD mi 

9Ca 


NO * NO «• 1 

909 


JO OF (NOp jKD ) • NO 

910 

139 

CONTINUE 

9 1 1 


theta c theta ♦ thetao 

9 i 2 

1 90 

continue 

913 


IF U «EU< 2 »0H« I «&Q« 

9 1 9 


NDP s NDP - NClft 

915 


N 1 = NOp 

9 1 6 


N2 s* NOp * NDV2 ♦ 2 

9 17 

C 


9 I 8 


NJ S N1 +| 

9 1 9 


N2 * N2 - 1 

9 20 


JD0F(N2,1) = -U00F1N1 t i | 

921 


UUOF ( N 2 ? 2 ) ■ JOOF ( N 1 ,2 ) 

922 

c 


923 


00 19l J 3 1 1 

929 


N i » N l + 1 

92b 


N 2 » N 2 * 1 

926 

c 


92 7 


JOOF ( N2 9 1 ) a ~JDQF(Nl ,1 2 

92a 


JOOF (N2 t 2) * J00F(Nl,2) 

929 


JOOF (N2 S 3) a «JD0F{Nl*3| 

9 3Q 

191 

CONTINUE 

93 l 


GO TO 1 l 29 

932 

t 122 

continue 

933 

C 


939 


nop 3 nop “ ncip 

93b 


N 1 3 NOP 

936 


N2 * NOP <► NDV2 ♦ 2 

937 

C 

- 

936 


N| • Nl ♦ | 

939 


N 2 B N 2 «* 1 

99Q 


JOOF<N2 9 li * -JOOFlNl,!) 

99] 


JOOF(N2 B 22 a -JD0F(N1»2| 

992 

C 


993 ’ 


00 1 I 23 J a 1 f NW| 

999 


N 1 85 N I ♦ 1 

995 


N 2 * N 2 - 1 

996 

c 


997 


JD0F|N2 f U 3 -JOOF ( N 1 e 1 > 

996 


JD0F(N2j2> * “ JOOF (N 1,2) 

999 


JOOF ( N2 | 3 ) • JOOF {Hi , 3 ) 

950 

1 123 

CONTINUE 

951 

C 


9 52 

1 1 2 9 

continue 

953 

C 


959 


N 1 a Nl * i 

955 


NDP = NOP ♦ nCIN 


NPI > GO TO H22 
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*456 
*457 
*450 
*459 
<460 
*46 \ 

*46 2 
*463 
*46*4 
*4 6 S 

*466 

*467 

*468 

*469 

*470 
*47 1 

*472 
*473 
H7H 
*47 b 
*476 
H 7 7 
*4 7 b 
*479 
*480 
<461 
*482 
*463 
*48*4 
•*48b 

*466 
*48 7 

*468 
*489 
*490 
*49 1 
*492 
*493 
*49*4 
*4 9 b 

*496 

*4 9 7 
*496 
*499 

500 

501 

502 

503 
50*4 
50b 
5C6 

507 

508 

509 

510 

511 

512 


120 CQNT 1 NUfc 
200 continue 
C 
C 

CALL WRITE (XYZ ?NOP 8 3»3HXYZ »KX» 

CALL WRITE ( £UL » NDP , 3 , 5HEULLK , K * > 

Call WRITE UYZ8,NDPb,3,*4HXYZB»KS) 

CALL WRITE UUt-a,N0P&,3,6HEULER6»KS) 

CALL WR 1 T 1 fl ( jDUF ,NDP , 6 ,*4NJD0F *KX ) 

C- BRING XYZ,XY2B and EUUiLULB tq6ETher also add the jdofb* 

c 

DO 2 1 0 I=*l»NUP0 
NDP A * nDP ♦ I 
DO 220 J»I »3 
XYZ(NDPa»J) » XYZBU ,J> 

220 EUL<nOPA»J) ■ EULBd.J) 

210 CONTINUE 

C- JDOF TAbLE (TOP CENTER POINT! 


C- 


225 


2^0 


NDP m NOP + 1 

NCNR = (NplNTS l 2 ) 1 *NC 1R ♦ HPJNTS U) 

NO b Ni) H 
JDOF (NDP * 1 } 8 NO 

all OTHER POINTS ON the bLADoER (ON THE MUST BAY* 
NCNR = NCNR - NCIR - J 

00 230 i®l »NC1R 
NCNR 83 NCNR ♦ l 
NOP * NQP t 1 

JD0F(NDp,21 a JOOF ( NCNR * 2 ) 

IF ( J *6T. N02 > 60 TO 230 
IF (ICE;**! «0R« 1*EO*N02> 60 TO 2 25 

NO » NO ♦ 1 • 

JDOF ( NO p , 1 | 3 NO 
NO « NO + | 

JDOF ( NDp * 3 J « NO 
60 TO 230 

continue 

IF (I • EQ • n NO a NO + 1 

IF (I «E 0 . n JOOF (NDP i 1 ) ■ NO 

IF (I .to* N( 42 > NO « ND + 1 

IF (1 « E Q * Nw 2 ) JD OF ( NOP » 3 1 ■ NO 

CONTINUE 

NCNR » NCNR + 1 - NCIR 
NDP s NDP - NC 1 R 
N 1 ■ N D P 

N2 * NDp + NO V 2 + 2 

N 1 « N 1 * l 

N 2 * N2 •=* 1 

JD0F(N2, l ) = -JD0F(N1 , 1 ) 

J00F(N2,2) = -J00F(N1,2) 

00231 I * 1 ,N« 1 
N 1 * N 1 ♦ 1 

N2 * N2 - 1 


C 


JOOF (N2 , n * - J 0 0 F ( N 1 * 1 ) 
J0GF(N2j2) = -JD0F(N1,2) 
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513 
5 l *1 

sis 

516 

J 517 
516 
51? 
52q 

52 l 

522 

523 
52M 
S2b 
526 

52; 

526 

52? 

530 

531 

532 

533 
53*4 

53b 

536 

537 
536 
53? 
S*40 
SHI 
542 
5*43 
5*H 


JOOF ( N2 , 3 ) b J00F(NX f 3) 

231 CONTINUE 
C 

N l a N 1 + | 

NDP » NOP ♦ mcIR 
C 

c* ALU OTHER BLaDOLR POINTS EXCEPT T HE LAST ONE. 
NKBAYl ■ nrbay - I 
DO 252 KK*2jNR*AY1 

NCNR ■ NCNR NP l NTS f g K + I ) • NCIR 
DO 2*40 1«1 *NCI h 
NCNR * NCNR ♦ 1 
NOP s NOP ♦ i 

JD0F(NDp,2) = JDOF 1 NCNR , 2 ) 

IF (1 *0T. Nw2 1 60 TO 2H0 

IF .OR* 1 « EQ . NQ2 ) 60 TO 2*45 

NO * ND ♦ J 

JDQF {NOp tl) s NO 

ND = ND «■ 1 

JDQF ( NOp » 3 ) a NO 

60 TO 240 
245 CONTINUE 

IF (I *EQ. U NO ■ NO ♦ 1 
IF (1 « EQ # I) JDQF { NOP * I ) « NO 
IF (2 * EQ « N02> ND = ND ♦ 1 
IF (I *£Q. N02l JOOF|NOP,3) » NO 
240 CONTINUE 

NCNR * NCNR - NCIR 
NOP * NOP - nCIR 
N 2 a ND'p 

N2 * NDp ♦ N0V2 + 2 
C 


54b 



N 1 b N l -fr 1 


546 



N2 e NZ - 1 


547 



UD0F(N2 p U ® 

- JOOF (N 1,1) 

546 



JOOF ( N2 , 2 ) ■ 

-JOOF ( N 1 ,2) 

54? 

C 




550 



00 2 *4 1 1 * 1 

,NW1 

552 



N 1 m N 1 + 1 


552 



N2 ■ N2 ** 1 


553 

C 




554 



J 00 FiN 2 »n a 

-JOOF (N l , 1 ) 

555 



JOOF ( N2 , 2 ) a 

-JDOF(Ni ,2) 

556 



JOOF < N2 » 3 ) b 

JOOF i N 1,3 ) 

557 


241 

CONTiNUE 


556 

c 



4 

55? 



NDP = NOP * 

NCIR 

560 


250 

CONTINUE 


561 

c 




562 

c- 


take care of 

THE LAST POINT WITH NO O.q.F* 

563 



NDP B NOP ♦ 

NCIR 

56*4 

c 




565 

c • 


FOR RGTaTjQNaL DO F AT THE BLADDER POINTS 

566 

c 




567 



NOP = NOP - 

NOPS 

566 



NDP b NOP ♦ 

l 

56? 



00 260 kK*1 , 

NRB A Y 1 



A2-33 


t>0 270 1*1 * N W2 

If <1 «E<*U N02> GO TO 26M 

IF ( l »NEc 1 ) GO TO 262 

NOP a NOP ♦ 1 

NO * NO ♦ 1 

J00F(NDp,6) a NO 

GO TO 270 
262 continue 

nop o nop * 1 

NO = NO ♦ 1 

JOOF(NOpiH) a NO 
NO * NO * 1 
JOOF(NDp»b) a NO 
NO * NO ♦ 1 
JOOF ( NOp t 6 ) a NO 
GO TO 270 
2 6 h continue: 

nop B NOP ♦ ] 

ND * ND + 1 
J00FlN0p,<O e NO 
N 0 ■ N 0 + 1 
J00FlNDp,5) = NO 

2 7 0 CONTlNUt 

NDP » NoP - NQ2 
Nl * NOP 

N2 * NDP NO V 2 -♦ 2 
C 

N 1 a N 1 + 1 
N2 a N 2 - I 

JOOFlN2, t ») « JOOF { N I t <♦ ) 

JuGF(N2,5) * JOOFINIjB) 
J00F(N2,6>»-U00F(Nl|6) 

C 

00 27i I » I ,NJ1 

N l a N 1 + I 

N2 * N2 - t 
C 

JD0F{N2 ( H) » JDQF ( N 1 , q) 

JOOF (N2 ,5 ) a JOOF (Nl 9 $ ) 

JOOF (N2,6) » -JOOF (Nl r 6> 

271 CONTINUE. 

C 

NOP s NuP + NCJ K 
260 CONTlNUt 
C 

c 

Call WRITE ( x Y2 ,NQPA . 3 ,3HXYZ ,KX > 

CALL WRITE (EUL »n0PA t 3 ,5HEULfcR »*X ) 
call writim < jduf , ndpa , 6 , mhjdof * xx > 
c 

c 

WRITE (NUTXYZI NDPA ,NCX , NOPa , nc J *NOPA I NCX 
rt'RI TE < NOTX YZ ) ( ( JDOFU # J J , I al ,NOPA ) , J»i ,NCJ» 

WRITE (NUTXYZ) (( XY2 u t J » , 1 *1 iNDPA) * jal ,NCX ) 
WRITE (NUTXYZ) (( EUU( I »ji I I »l sNDPA) ,J»1 ,NCX) 
C 



62 ? 

628 

62V 

63$ 

631 

632 

633 
63H 

635 

636 

637 
636 
63? 
6HQ 

6H1 
6 Hi 

6R3 

6RR 

6^5 

6H6 

6R7 

6R8 

6R9 

650 

651 

652 

653 
65*4 
65b 

656 

657 
6 5 8 
65y 
660 
66 1 
662 
663 
66R 
66b 
666 
667 
66b 

669 

670 

671 

672 

673 
67*4 

67b 

676 

677 
67b 
679 
6 80 
681 
682 
683 


C- ELEMENT* IN thE FIRST BAT 

E PENTaS 


N J EE 1 

NCOMfrlT * NELEMi ( n - 
DO 3 1 0 1*1 » N C 0 D N T 
NJ e NJ + NCIR * 1 

N 1 * N J 
N2 s NO «■ 1 
N 3 * N J 

NR ■ NO - NCIR “ I 


32Q 

310 


N 5 * NO - N C i R 
Nfc * NO « NC l R 

- 1 

DO 320 O^IeNCIRI 

NEL * NEL + 
N 2 * N 2 ♦ 1 

N 3 ® N 3 ♦ 1 

i 


N 5 s N 5 + 1 
N 6 s N 6 + 1 



ICONF tNEL 8 l ) 

= 

Nl 

ICONF < NEL ,2 ) 

B 

N2 

ICONF ( NEL ,3 ) 

a 

N 3 

ICONF t NEL , R ) 

s 

NR 

ICONF CNfcL 5 5 ) 

at 

N5 

I CONF tNEL t 6 J 

® 

N 6 

continue 

CONTINUE 



TETRA t T OP ) 



NJ = NO + N C 1 R 

* I 

NR * NO 



Nl ■ NO - NCIR 

« 1 

N 2 * N 1 

NB J a NTPF + 

l 


NB 2 ■ Nbl 



DO 330 I*l»NclRI 

NEL * NEL + 
N2 * N2 + t 

1 



l 


N3 « N2 ♦ 1 
ICONF <NEL & 1 ) * N 1 
ICONF tNEL t 2 ) = N2 
1 CONF <NEL ,3 | s N3 
ICONF (NeL »R ) » N H 
NELS ■ NELS + i 


icons <n£ls» i j * nr 

ICONS (N£LS«2 ) s N2 
1 CONS ( NtLS * 3 ) * N 3 
NEL.B = NEL8 + 1 


NB2 *, Nb2 ♦ l 
NB3 * Nb2 + l 
l CONb 1 NEL8 si) * NB1 
iCONBCNELB # 2 ) 38 NB 2 
l CONS CNlLB * 3 1 ■ NB3 
330 CONTINUE 

N J <s NO <= 1 N£ I R ♦ | ) 

NO * NO « (NCJR^I^ * NCOUNT 
NO 1 o Nj 

N 02 * NPINTS <l> * IN6INTS i2)*NCJR) 
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68M 

685 

66 * 

68 7 
6 fa B 
669 
690 

69 l 

692 

693 
69M 
6 9 & 


C- K£A0 IN ThE bQTTQM AND TOP fcLEMtNTS IDENTIFICATIONS* 
head <&»30Q2i ecmiqi ,elmid2 
C- SUTTOM ELEMENT of the second bay* 
head (5,3002) ELMIDI .ELHI02 
c 
c 

IF IeLMidi ,n£* SHPENTA) go to 3H0 
N l a N 0 I 

(42 = NO l + NClH ♦ i 
N3 o NJ2 

DO 3 faO IMiNCIM 
NEC b NEL ♦ i 


696 

697 
696 
699 


N i * N 1 ♦ 

N 2 = N2 * 
N 3 ® N 3 ♦ 
NM » N I ♦ 


70C 

701 


NS e N2 ♦ 
N6 « N 3 ♦ 


702 


ICON? (NtL» 1 ) « N| 

703 


ICONF (NEL, 2 ) b N 3 

70 M 


IC 0 NF{NEL, 3 ) 0 N 2 

70 S 


1 CONF < NEL 1 ■ NM 

706 


ICONF <NEL , 5 ) » N 6 

707 


I CQNF < NEL , 6 ) a NS 

706 

3&0 

continue 

709 


NJJ a NS 

710 


NU 2 ■ N 6 

7 1 l 


GO TO 370 

7 l 2 

3 M 0 

CONT 1 NUE 

7 I 3 


N I s NUi + NClH + 

7 l M 


N 3 3 N J 2 ♦ NCJH ♦ 

7 15 


N 5 a NUi 

7 16 


N 7 s N 02 + l 

7 1 7 


DO 360 1 = 1 »NC I Hi 

7 18 


NEL = N£L + 1 

7 19 


N I - N 1 * 1 

720 


H 2 « N 1 + l 

721 


N 3 3 N 3 ♦ l 

722 


NM * N 3 - I 

723 


NS * NS ♦ 1 

72 M 


N 6 s N 5 + 1 

7 2 S 


N 7 ® N 7 ♦ l 

726 


NS = N 7 - I 

727 


ICUNF<NEL,l ) * NJ 

728 


ICONF <NEL # 2 ) a N 2 

7 29 


IC 0 NF(NEL, 3 J a N 3 

730 


ICONF { NtL J B NM 

731 


I CONF (NEL ,5 ) ■ NS 

732 


I C ON F ( N E L , 6 ) 0 N 6 

733 


ICONF {Nel ,71 a N 7 

73 H 


ICONF ( N E L 1 & ) * N 8 

7 3 S 

360 

CONTINUE 

736 


NJI s N 2 « NC I K 

737 


N J 2 * N j - NC I H 


< 


73a 370 CONTINUE 

7 39 C- In BETWEEN EeEHENTS IN THE SECOND BAY {EXCLUDING TOP ♦ BOTTOM ) 

7 M 0 C 
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7 H 1 
7H2 
7 M 3 
7Hh 
7MS 
7^6 
7 M 7 
7MB 
7 M V 

7 bO 
7S1 
7&2 
7 5 3 
7b M 
7 b 5 
7 S 6 

757 

758 

75? 

76Q 

761 

762 

763 
7 6S 

76 S 

766 

767 
7 6 8 
76V 
7 7 C 

771 

772 

773 
7 7M 

77b 

776 

171 
778 
77? 
780 
78 1 
782 
76 3 
7 8M 


NcOUNJ * NfcLEM l < 2 I « 2 
NI « NJi ♦ NtI K ♦ i 
N 3 = NU2 * NCl^ * 1 
NS = N J | 

N7 b NJ2 + | 

00 380 I*1*NC0UNT 
00 3?0 0=1 »NC I R i 
NEt S NtL ^ 

N| a Nl + | 

N 2 « N I * i 
N3 * N3 + 1 
NM = N3 « | 

N5 « N5 + l 

N6 = Nb 4> 1 

N7 » N7 + l 

N 8 * N 7 — I 

ICONF <NEL,'I » » Nl 
ICONF ( N£L » 2 I a N2 
ICONF i N£L j 3 ) o N 3 
ICONF IN£L ,M ) a NM 
ICONF ( NtL > 5 I a NS 
1X0 NF(NEL v 6) a N 6 
I C ONF I NtL » 7 ) B N 7 
I CONF <NtL,8 ) * N6 
370 CONTINUE 

Ni B N| - NCIRJ ♦ NCIR ♦ 1 

N 3 = N3 - NC I R i «• NCIR 

NS a NS - NCiRl ♦ NCIR ♦ I 

N 7 b N 7 ™ NCiRl + NCIR 

Njje;N2™NCl^ 

N J2 a N3 - NC1« “ I 
3SQ CONTINUE 
C 

C- FOR THE LAST ELEMENT 

NKOUnT * NPINTSI2) 

IF (ELMjOi *EQ* 5HPENTA) NKOUNT * NKOUNT 
IF INKOUNT .G*. NPIntS | 3 ) > NJ3 a NUl ♦ 

IF (NKOUNT *LT 6 NPlfcjS (3)) NU3 a NJ2 + 

IF (NKOUNT ,E«» NPJNTS (3)) GOTO MOO 

C 

Nfll a NTPF + I ♦ i 
N83 o NTPF ♦ NCIR ♦ I 
N l a N J I 
NZ a NU2 


76b 

786 

767 

786 

78? 

7?0 

7? J 
7 ?2 
793 
7 9 M 

795 

796 

797 


N 3 = N J 3 

00 MjO I»| e N C I R 1 
NEL = NfcL + l 
N| e Nl ♦ 1 
N 2 * N 2 * 1 
N 3 a N 3 * I 
NM * Nl + J 
NS a N 2 «■ l 

N6 at N3 * 1 

ICONF (NtL „ I ) a N 1 
I CONF ( NtL b 2 ) s N 3 
IC0NF(NEL»3l b N2 
iCONFtNELjM) « NM 


- 1 

NCIR + I 
NCIR 



79# 

799 

800 
*801 

802 

603 

8Q4 

806 

806 

607 

808 

809 

810 
8 i i 
8 1 2 

813 

814 
61b 

816 

817 

618 

819 

820 
821 
822 

823 

824 
8 25 
826 
827 
82b 

829 

830 

831 
332 
833 
83*4 
83b 

836 

837 

838 

839 
84Q 
8H 1 
8 H 2 

843 

894 

84& 

696 

697 

898 

899 

8 &0 
651 

8b2 

853 

85<4 


1C0NF ( N £' L | 5 1 » N 6 
1C0NF (NtLj6) « NS 

IF (NKQuNT ,LT* NPInTS,( 31) 60 TO 395 
MS! a N 3 

NS2 * N2 

N S 3 * N 5 

NS9 * M6 

NB1 • Nbl ♦ l 

NB2 a Nbl - 1 

N# 3 * N#3 ♦ 1 

NB9 * N63 ♦ 1 

60 TO 397 

395 continue 
MSI s N 1 

NS 2 a N 3 

NS3 v N 6 

NS9 a N9 

N B 1 = N £J 1 ♦ \ 

N02 a Nbl - 1 

NB3 a Nb 3 ♦ l 

N89 a N# 3 + i 

397 CONTINUE 


nelb * NELB ♦ 

1 


ICONS INELB t l ) 

» 

NBI 

1 CONb ( NEE8 t 2 ) 

C 

NB 2 

1 C0N5 ( NELB » 3 ) 

fit 

NB3 

ICON# JNEE8»9) 

m, 

N84 

nee 5 ° nels ♦ 

i 


icons (nees » i > 

a 

NS 1 

ICONS<NEES >2) 

m 

NS2 

I CONS ( NELS » 3 ) 

S3 

NS3 

IC0N5<NLLS » 9 ) 

u 

NS4 


9 i o continue 

IF (NKOuNT ,GT» NPINTS (3)> N02 a N02 
IF (nKOuNT , L T • NPInTS (31) N02 a Nu3 
60 TO 930 
9CC CONTINUE 

N03 * N J2 + NCJR 
NU9 = NO I ♦ NC I K ♦ 1 
N 1 ■ NO 9 

N9 * NU3 
N5 = NO 1 

N 8 « N 02 

UO 920 I = | tNClKI 
NEE ■ N t L + 1 
N l — N 1 + 1 

N2 a N1 * 1 

N9 = N9 ♦ I 

N 3 a N 9 + 1 

Nb = N5 ♦ 1 
N 6 = N 5 ♦ 1 
N 8 a N 8 ♦ 1 

N 7 “ N 8 ♦ 1 

1C0NF (NtL , 1 ) a N 1 
1 CONF ( NEL b 2 ) a N2 
1C0NF(NEL,3| a N3 



A2-38 


6 5$ 

ICONF <NEt ,9 ) 

« NR 

856 

ICGNFWL, 5) 

a NS 

657 

ICONF <NtL ,6 ) 

s N 6 

650 

ICONF < NtL , 7 ) 

B N 7 

859 

ICONF (NLL ,8 ) 

a N 6 

860 

N 5 1 * N2 


66 1 

NS2 « N 1 


662 

NS 3 » N9 


863 

NS9 « Ni 


869 

NELS “ NELS + 

I 

665 

ICONS (NELS • I ) 

B J 

866 

I CONS ( NELS t 2 ) 

■ N52 

867 

ICONS (NELS » 3) 

a NS3 

868 

ICONStNELSjR). 

a NS9 

B69 

NELB ■ NELB ♦ 

I 

8 7 0 

NBl a Nbl ' ♦ 1 


67 1 

Nb2 « NB I ♦ 1 


672 

NB3 « NB3 ♦ 1 


873 

NB9 a NB 3 + 1 


67 9 

ICONB i NEL6 » 1 ) 

» NB 1 

875 

l CONb ( NELB p 2 ) 

« NB2 

876 

(CQNMNELB s 3 ) 

« NB3 

877 

ICONB^NEUB »9) 

a NBR 

678 

920 continue 


879 

N J 2 a N 0 3 


880 

NJ3 a Njq 



86 X 
882 
883 
ass 

88b 

866 


930 CONTINUE. 

C 

C- ALL OTHER ELEMENTS FROM THIRD TO LAST BUT ONE HAT* 
C 

' NtJTl * NTPF * I + I 
N6T3 * NTPF + NClR * | 


687 

688 


N81 a NBTI 
NS3 a N6T3 


869 

69C 


NJI ■ NU3 ♦ NC iw ♦ 1 
NJ2 » NU2 nCIR 


891 

892 

8 9 3 
899 
898 
896 


DO 5Q0 J I «3 i n RB A Y | 

NBJ * N8I ♦ nCIR 
NB3 « NB3 ♦ NCIR 
READ ( 5,3002) ELM I 01 ,£LMID2 
C- FIRST (BOTTOM jELEMENT OF EACH 

IF IELMJOI'.nE* 5 H P E N T A ) 60 TO 510 


697 
8 9 6 


NUU1 a NJJ 
N J J 2 s NJ2 


899 

NJU3 ■ 

NO I 

♦ NC I R 

900 

N l 

a N J J 1 


90 1 

N 2 

» NUU2 


902 

N 3 

b N J J 3 


903 

DO 

520 

i«l 

*NCl H l 

909 

nel 

a NEL 

♦ I 

905 

N l 

• Nl 

+ 

l 

906 

N 2 

a N 2 

♦ 

I 

907 

N 3 

a N3 


l 

908 

Nr 

B N 1 

♦ 

1 

909 

N5 

a N 2 

♦ 

l 

910 

N6 

a N3 

♦ 

I 

9 I 1 

iconf (nel i 

1 > a Nl 




9 U 

913 

914 

915 

916 
91? 

91b 
9 19 

920 

921 

922 

923 

924 
92b 

926 

927 
92b 

929 

930 

931 

932 

933 

934 
93b 

936 

937 
93b 

939 

940 
94 1 

942 

943 

944 
94b 

946 

947 
94b 

949 

950 

951 

952 

953 

954 

955 

956 

957 
95b 

959 

960 
96 1 

962 

963 

964 

965 

966 

967 
96B 


1C0NMN£L,2> a N3 
IC0NFiN£L»3) = N2 
ICONF 1 *4 1 L * 4 ) a N4 
ICONF < N£L } & ) * N 6 

ICONF<N£l_ ,6) * Mb 
b20 CONTINUE 

NJCi ■ NJJ3 
NJC2 ■ NJJ2 
(jO 10 540 

sic continue 

NJJJ = NJl 
NJJ2 * N 0 2 
N J J3 » nJ 2 + NCIR 
NJJ4 ■ NJl + NCIR 
Ni a N J j 4 
N 4 « N J j 3 

N5 = N Jo 1 
N b b N J J 2 

DO 530 .1 M »Nc 
NEL * N£L ♦ 1 

N l a N 1 <* l 

N2 « Nl ♦ | 

W 4 « N 4 ♦ 1 
N3 « N4 * | 

NS s NS + 1 
N 6 » NS ♦ 1 
N Q e N 8 * 1 
N 7 a N 8 * 1 

ICONF (NEC , 1 ) b Nl 
lC0NFtNEL,2) = N2 
IC0NF<NEL,3) a N3 
ICONF (NeL,4) s N 4 
ICONF < NEE ,5 ) a Nb 
lCONFlNEL.,6) b N 6 
ICONF tN EL *7) a N7 
IC0NF<NEL # 8) a Ne 
S 30 CONTINUE 

NJCI 3 NJJ4 
NJC2 * N J J 3 
S4C CONTINUE 
C 

IF ( NELEM 1 I J I ) «CE* 21 00 TO JD&O 

c- in betacEn elements elements of each 

NCOUnT * NELEM Mil) - 2 
00 550 1=1 1 NCOUNT 
N JC 3 3 NJC2 + NCIR 
N JC 4 * NJCI + NCIR 
NJ = NJc4 
N 4 a N JC 3 
N b = N J C 1 
N B = N J C 2 

DU 560 J«1 1 NC I R l 
NEL a NEL ♦ l 

N J a N 1 «■ 1 

N 2 ** N 1 + l 
N4 a N 4 ♦ I 
NJ a N4 + 1 
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9 6 y 

97o 
9 7 i 

972 

973 
97s 
97b 

976 

977 

978 

97 V 

98Q 

981 

982 

983 
98i* 

9 8S 

986 

987 
980 

989 

990 

991 

992 

993 

994 

9 VS 
996 

99? 

998 

999 
1 000 
1 001 
1002 
1 003 

loos 

1 QQS 
1 006 

1007 
1 QOS 

1009 

1010 

1011 

1012 
1013 
101S 
10 18 
1 Q 1 6 
1017 
10 18 

1019 

1020 
1021 
1 022 

l 023 
102S 
1 Q2S 


1*5 a N5 ♦ 1 

N6 ® NS ♦ 1 
N 8 a N 8 «■ 1 
N7 * N6 + 1 

JCQNF IN£L * 1 ) ** Nl 
1C0NF INfcL 9 2 ) « N2 
1C0NF(N£L # 3) » N3 

ICONFtNtL 8 S) m NS 
l CQNF 1 NtL y & ) a N5 
ICONF < NlL ,6 ) b N6 
lC0NF<NtL,7> o N 7 
I C QNF i N£ L » 8 ) * N8 
S6C continue 

NJCi ** much 

N JC 2 B NJCS 

S&0 continue 

c 

1050 CONTINUE 

IF 1 NELEMUlIl aUTfl 2) t»0 TO &8C 
C- FOR THE LAST ELEMENT IN EACH B AT * 

L 

NKOUNT - NPINTS(ZI) 

IF tEOrtm » EQ • SHPENTA) NKOUNT « NKOUNT - 
IF (NKOUNT ,gT. NPInTS m-M>» N JC 3 » NJCI 

IF (NKOUNT * UT • NPINTS Ul + 111 N JC 3 ■ NJC2 

IF (NKOUNT ,eQ* NPINTS til-Mll 00 TO &6& 

N 1 a NJCI 
N2 * N JC2 

N 3 « N J L 3 

00 5?0 IMtNCl^l 

NEL = NfcL * l 
N J = N I ♦ 1 

N2 ® N2 ♦ 1 

N3 ss N3 + 1 

N S « N I * 1 

N& « N2 + 1 

N 6 a N 3 ♦ I 

ICONF 1NEL» i ) » N 1 
lC0NFtNEL»2l b N3 
I CQNF l NEL 8 3 ) = N2 

IC0NFINEL 9 S1 a NS 
iCONF(N£L,b) « N6 

1 CO NF ( NEL | 6 ) a NS 

IF (NKOuNT * L T • NPInIS (II-MJl 00 TO bS5 

NS 1 a N3 

NS2 a N2 

NS 3 = NS 

NS S * N6 

00 TO &S7 

&SS continue 

N51 s I*} 

NS2 a N 3 < 

NS 3 = N 6 
N5S * NS 
S5'7 CONTINUE 

nels s NELS «• l 
ICONS < NELS t 1) a NS 1 


l 

NC I ft 
NCIR 



026 

027 

02a 
02 <? 
030 
0*1 
0*2 
0*3 
03^ 
0*5 
0*6 
0*7 
038 
0*9 
Q4Q 
OH 1 
0H2 
043 
0 4 4 
Q45 

046 

047 

048 

049 

0^0 

051 

052 

053 
059 
05b 
0 5 6 

057 

058 

059 

060 
06 i 
062 
06 3 
064 
G 6 b 
066 

06 7 

068 
G 6 9 
070 

07 1 

072 

073 

074 
C 7 S 

076 

077 
076 
D 7 9 
060 
081' 
082 


ICONS <NLLS >2 ) ■ NS2 
I CONS t NCL S * 3 ) « NS3 

lCQNS<NfcLS,4) * NS4 

HE LB ® NELB ♦ i 
NB 1 s NBl + 2 

N82 * Nbl • 1 

NB3 ■ Nb3 + l 

NB4 a N63 ♦ 1 

1 CONB (NfcLB • 1 ) ® NBl 
1 CONb ( N£LB » 2 > ® NB2 
ICONS CNELB »3) ® N&3 
IC0NB<N£L8»4J ® NS4 
& 7 0 CONTINUE 

IF (nKOUNT ,GT. NPINTS IliMM NjCi * NJC3 
IF (N^OgNT ,lT* NPJNTS (!!♦!)) NJC2 * NJC3 

GO TO 580 
56b continue 


N JC 3 « NJC2 + 

NC 1 R 

NJC4 ■ NUC l + 

N C 1 R 

Ni a NJC4 


N 4 “ N si C 3 


NS « N J C 1 


N 8 « N JC2 


00 S90 !■ I ,NC I R l 

NEL * N EL + 1 


N I EB N i ♦ 1 


N 2 * N 1 * 1 


N 4 ® N 4 ♦ 1 


N3 » NH ♦ 1 


NS » NS + I 


N 6 s NS * t 


Nti ® N 8 + t 


N 7 a N 6 + I 


ICONF<NEL f l) = 

NI 

IC0NMNEL,2) a 

N2 

ICONFINEL,*) = 

N3 

I CONF ^ NEL p 4 ) B 

N 4 

1C0NF < N E L # 5 ) « 

NS 

I CONF t NEL » 6 ) * 

N 6 

IC0NFCNEL.7) b 

N7 

IC0NFlN£L # 8) * 

N 8 

nels ® NEL 5 + 

1 

MSI a Ni 


NS2 ® N 4 


NS* a N4 «■ 1 


NS4 ® Ni ♦ l 


I CONS ( NELS » 1 ) 

a NS 1 

ICONS «NELS t2) 

® NS2 

1 cons i nels » 3 1 

a NS3 

IC0NS<N£LS ®4 ) 

a NS4 

NELB * nELB «• 

1 

NB1 e NtH ♦ i 


NQ2 * Nbl - l 


N & 3 « Nfi 3 ♦ 1 


NB4 a NB3 ♦ 1 


l CONB l NLL6 t 1 ) 

* NB 1 

I CONB i NELB 9 2 > 

* NB2 
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083 

084 
08b 
086 

087 

088 

08V 

OVQ 

ovi 

0 V 2 

093 

Q9R 

098 

096 

097 

098 
09V 
100 

101 

102 

103 

104 

10b 

106 

107 

1 08 
1GV 
1 10 
l 1 1 
1 12 

113 

114 
lib 
116 
117 
1 18. 
1 IV 
120 
121 
122 

123 

124 

12b 
■ 1 26 

127 

128 
1 2 V 

130 

131 

132 

133 
13S 
1 3b 
13b 

137 

138 
13V 


1C0NB (NfcLB ,3 ) * NB3 
lCGN6<Nfc.LB»4) *= N84 
590 CGNT 1 NU t 

NJCl ■ N JC 4 
N3C2 a NJC3 
58 Q CONTINUE 

Nfll e Nbl - nCIRI 

N83 = Nfa3 " N C 1 R 1 

NJ1 s N J 1 <9“ nC1K*NP1nTS (ID 

NJ2 « N J 2 * NCiR*NPlN|5 ( l 1 ♦ 1 ) 

500 CONTINUE 

O FOR JH£ OAST <QNE ELEMENT) BAY • 
NJ3 = Nij 2«* Ncl« 

N 1 s» N J 1 
N 2 « N32 
N3 ■ NJ3 

NB1 h N61 ♦ NC Jft 
Nb3 « Ntt3 ♦ NCIR 

uo 6 } 0 1*1 » N C 1 N 1 

NEL * NtL ♦ i 
N l « N 1 ♦ 1 

N2 » N2 + 1 

N3 ■ N3 > 1 

NR « N 1 + l 

N5 « N 2 + 1 
N 6 s N 3 + 1 


iCONplNELpll « 

N 1 

IC0NF<NEL,2) » 

N 3 

lC0NF(NELp3) ■ 

N2 

1 C ONf C N£L * 4 ) * 

NR 

ICONFtNELpb) b 

N6 

I C ONF 1 NEL , 6 ) b 

N5 

«EL 5 ® nELS + 

1 

NS 1 = N3 


NS2 - N2 


ns 3 s Nb 


N S 4 = N6 


ICONS CNELS » 1 1 

* NS 1 

1 CONS 1 NcL5 »2 ) 

- NS2 

ICONS < NlLS * 3 ) 

* N S 3 

lC0N5lNtLSt4) 

= nsr 

NEC 8 * NEL 8 ♦ 

l 

N8i « NB1 ♦ 1 


N 8 2 = N 8 1 - 1 


N 8 3 ■ Ny 3 + l 


Nb4 s N8 3 ♦ l 


ICONB i NELB t 1 ) 

* N0 1 

1 CONS 1 NEL 8 » 2 1 

® N&2 

ICONS! N t L8 » 3 } 

= N 8 3 

ICONB ( NELb 1 4 ) 

* NS 4 


6 1 0 CONTINUE 

CALL wRMIM ( j CONF p NEL , 8 > HH JF LO i l 500 I 

c----- — - — — " 

c 

C- «KITE FlUIO PROPERTY ON TAPE (INTERIOR) 

WRITE (NUTfcL *30011 A N 6 ( 2 ) 

WRITE tNUTEL, 30021 AN* l 7 1 « ANb t 9 | » ANb l 4 ) , ANb < 1 l 1 » AN6 t 4 1 
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l 140 
l 141 
1 t 4 £ 
l 143 
1 144 
1 14& 
1146 
» I 4 7 
1 146 
1 149 
1 I BO 
1 1B1 
1152 
1 l b 3 
1 154 
1 155 
1 156 
l 157 

1 158 
l 159 
U60 
1 161 
l 162 
1 163 

l 164 
1 1 65 
1 166 
1167 
l 166 
1 169 
1170 
U71 
l 172 
I 173 
1 174 
l 1 7 5 
1176 
I 177 
1 1 7 b 
1179 
1160 
1 1 SI 
1 1 B 2 
1 163 
1164 
1 165 
1 166 
1 187 
1166 
I 1 69 
1 1 9 C 
i 1 91 
l 192 
I 1 93 


WRITE ( NUTEL , 3^03 t AN& I 1 ) » RHOF » ANS < 2 1 , BLKM 
c- WRITE fluid el £ «e nt s on TAPE 
DO 4Q5 L 3 l»NtL 

405 WRITE <NUTEL,30Q6j L, llC0NF(LiUl*U B l»6l 
wfRJfE ( NUTEL « 3306 1 UNT6IU1 fU al » 2 > 


c- write property of gravity element on Tape 
c 

CaEL A R l T I M ( ICONS »NEL5*4,4HjSUR t KS) 

WRITE (NUTEL,300ll An6(31 

WRITE InUTEL»3002J AN 6 ( 4 1 » A Nb ( 9 > » AN6 ( 4 j , AN 6 ( 4 1 » ANb t 4 > 
WRITE (NUTEL *30031 AN£ ( 1 ) * RHOF 

WRITE 1 NUTEl *3 003 1 AnS < 3 1 t G VA . AN5 ( 4 ) ,G V Y *An5 ( 5 1 *GVZ 
0- WRITE GRAVITY ELEMENTS On TAPE 
DO 7 1 0 Li = i„NELS 
L « l! + nEl 

7 IO WRITE IWUTEL*3p06l L, ( I CONS < L 1 » U J » J« 1 , 4 1 
WRITE (nUTEL» 3006) (INTG(Ul»U a l*21 


C 

CCCCCC- — — • — — — — — — — — 

C- WRITE PROPERTY OF TRIANGULAR BLADDER ELEMENT (NCiRl) ON TAPE 


CALL WR1TIM ( ICONB ,NELB *4 ,4HJBLU*KS 1 
WRITE (NUTEL *30011 AN615) 

WRITE {NUTEL,3002l An 6 ( 8 1 , AN6 l * ) , AN 6 l 4 ) , AN6 < 4 1 » AN6 1 4 1 
WRITE (NUTtL, 30031 ANS ( U , RROBLD * ANS < 6 ) *EBLD , AN5 ( 7 1 , ANUfaLD 
WRITE (NUTfcL, 30031 AN§(81 ,TBLD*AN5t9) ,TBLD,AN5U01 fTBLO 

c- write bladder triangular element on tape 

do 725 LI - i * NC I R 1 
L » L 1 + NEL ♦ NEL5 

725 WRITE iNUTEL, 30041 L, (ICONB ( L 1 , 3 1 . J« 1 , 3 1 , T BLD , T bLD . TBLD 
WRITE (NUTEL,3006) U NT G ( J 1 * J' l t 2 1 
IF (nEL 0 .EQ, NCIR1) 00 TO 727 
C 

C - WRITE BLADDER uuad. property on tape 

WRITE ( uuTEl * 300 1 1 A N 6 ( 6 1 

WRITE ( N U T E L * 3 0 0 2 1 AN6181 » AN6 ( 9 1 * AN6 ( 4 } * AN6 ( 4 ) * AN6 ( 4 I 
WRITE (NUTEL,3CQ3) ANSI 1.) ,RHOBLD ,AN5(6| # E»L0 # AN5(7) .ANUfiLu 
write <NUTEL,3003) ANS < 8 ) * T BLD * AN5 < 9 1 » T 8L0 » AN5 U 01 » ToLD 
C- write bladder rUaD. elements ON TAPE 
, DO 726 Ll * NClK,NELb 

L a L 1 + NEL ♦ NEL5 

72b WRITE (nUTEL. 30051 L, (I CONb ( L 1 » J ) » J* 1 , 4 1 , TBL U .TBlD , TBlO 
WRITE l NUTEL * 3006 1 I INTL ( J 1 * J* 1 *2 1 
111 CONTINUE 

WRITE (NUTEL , 3.C0 l I AN61U 
C 

C--- 

c 


return 

end 
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Ap pendix - B1 


EXPLANATION OF 

INPUT TO STATIC FREE SURFACE PROGRAM 


Card 
Nos . 

Input 

Explanation 

1 

Run no. j cols. 1-6; name, cols. 11-28 

Three cards to satisfy 

2 

Title 1, cols. 1-78 

subroutine "START" 

3 

Title 2, cols. 1-78 


4 

Order of the polynomial, acceleration due to 
gravity (in/sec^) , ullage pressure (psi) , 

Format (15,3017.8) 


ratio of volume fill to total volume 


5 

'STOP 1 , cols. 1-4 

End of data 
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i 

? — 
3--- 

A 

S--.-- 


bF.PbbU S I NGi"? 

- FRtt bUKp ACt F0« ^LKCUHY Si osh---------.- — PtnCtAif FULL 

G K A V I f y F 1 E L 0 - -•> - 0 , 0 0 3 >J 6 

* O*003«fc 0*10 U.60 

— STOP 






5? «4 
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Ap pendix - B2 


EXPLANATION OF 

INPUT FOR' DYNAMIC ANALYSIS PROGRAM 

, i 


Card 
Nos . 

Input 

Explanation 

1-31 


Control cards (if needed) 

32 

33 

34 

Run no., cols. 1-6; name, cols. 11-28] 

Title 1, cols. 1-78 r Three cards to satisfy subroutine 

Title 2, cols. 1-78 J "START" 

35 

Users comment cards 

Any no. of cards, the last card 
must be zeros. Cols. 1-10. 
Subroutine "COMENT" 

36 

1 INITIL 1 or 'NOINIT 1 , cols. 1-6 

To initialize or not to initialize 
the reserve tape 

37 

* GNXYZ 1 , cols. 1-5 

To call subroutine "GNXYZ3" 

38 

‘XYZEUL* , cols. 1-6 

To call subroutine "XYZEU3" 

39 

Radius of the sphere, cols. 1-10 

Format (F10.0) 

40 

2 x mass density, 2 x bulk modulus , 
of the fluid and gravity components 
in x, y, and z directions 

J 

Format (5F10.0) 

• 

41 

2 x mass density, 2xYoung's modulus,. 
Poisson*s ratio, and thickness of 
the bladder 

Format (4F10.0) 

• 


42 Number of circumferential points and Format (215) 
number of radial bays in the 
half space model 


43 

44 

45 

46 


Name, no. of rows, no. of cols, of 
the matrix 

The coefficients of the polynomial 
assumed for the free surface 
Zeros 


Subroutine "READ" 



Name , no . of rows , 
the matrix 

bay radii 

Zeros 


no . 


of cols, of 


Subroutine "READ" 

. 


Gig 

* 
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52 

53 

54 
58-70 

71-81 

82 

83 

84 

85 

86 

87 

88 

89 

90 

91 

92 

93 

94 

95 

96 

97 

98 

99 


Name, no. of rows, no. of cols, of 
the matrix 

Number of points on each vertical 

line 

Zeros 

X-coord. of the points on each 
vertical line from inside-out 

Identification of the first and the 
last element in each bay from bottom 


to top 



'FINEL', 

cols . 

1-5 

'MODES' , 

cols . 

1-5 


Name, no. of rows, no. of cols, of 

the matrix 

Zeros 

Number of modes wanted 

Number of modes used 

2 

Shift value for 0 ) (convergence wi 
be about this value) 

No. of maximum iteration allowed 

'PUNCH* or 'NOPNCH* , cols. 1-6 

'MECHEQ', cols. 1-6 


Name, logical tape no., run no. 


Last joint with only 3d.o.f., last 
d.o.f. number 

Name, no. of rows, no. of cols, in 

the matrix 

Zeros 

Name, no. of rows, no. of cols, in 
the matrix 

Column no. of rigid body matrix 
(non- zero) 

Zeros 

Nans, logical tape no., run no. 


V Subroutine -"REA DIM" 


Format (8F10.0) 


Format (2(A6,4X)) 


Calls subroutine "FINEL" 
Calls subroutine "0YM0DE" 


Matrix of assumed mode shapes (if 
any) . Subroutine "READ" 

Format (10X,I5) 

Format (10X,I5) 


Format (10X,E17.0) 

Format (10X,I5) 

Option for punch output 

Calls subroutine “MECHEQ" to calcu- 
late mechanical equivalent 

To read x, y, z locations of the 
points. Subroutine “READ" 

Format (215) 


U^ocation of the reference point. 
^Subroutine "READ" 

“'The elements of the matrix indicate 
which cols, of the rigid body 
‘matrix are non- zeros. 

Subroutine "READIM" 

Degrees of freedom matrix. ' 

. Subroutine "READ" 
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100 

Name, logical tape 

no . , run no . 

/ 

Euler angle matrix 
Subroutine "READ" 

101 

Name, logical tape 

no . , run no . 

Mass matrix 
Subroutine "READ" 

102 

Name, logical tape 

no . , run no. 

Modes matrix • 
Subroutine "READ" 

103 

* PLOT ' , cols. 1-4 



Option to plot the mode shapes 
(NO card, no plot) 

104 

105 

106 

Name, no. of rows, no. of 
the matrix 

Node no. to be '•plotted 
Zeros 

* 

cols, in 

■+ 

Joint numbers which are to 
i be plotted 

[ Subroutine "SEADIM" 

107 

Name, no. of rows, 
the matrix 

no. of 

cols, of 

x, y, z location of the points 
Subroutine "READ” 

108 

Name, no, of rows, 
the matrix 

no. of 

cols, of 

Euler angles of the points 
Subroutine "READ" 

109 

• Name, no. of rows, 
the matrix 

no . of 

cols, of 

Frequencies of the system 
Subroutine fr READ fr 

110 

Name, no. of rows, 
the matrix 

no. of 

cols, of 

Modes of the system 
Subroutine "READ" 

111 

Name, no. of rows, 
the matrix 

no . of 

cols, of 

d.o.f. matrix of the system 
Subroutine "READ" 

112 

'START*, cols. 1-5 



Preparatory to end 

113 

'STOP', cols. 1-4 



End of data 


114-118 


Control cards 
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1 

2 
3 
A 
5 


- R U N , 0 M08GRI ,HMo'1HEA8 704, S-BULTMOOOOO, 30, 2000/ 1200 BULTMAiM BIN S-190 


6 

=US*ER. 

-BOX PC H , S -195 -BUI. T • 

? 

= FREE 

T P F $ . 

3 

=ASG,T 

TPF$, F2/1/PDS/500 

9 

=ASG,t 

^,11,17701 

to 

= C 0 P Y , C 

> r. ,tp-$. 

It 

= F R F F 

T. 

12 

= A SG * T 

i . , f 4 o / 1 /posy io 

13 

= A S G , T 

2 . , F40/1 /POS/ 10 

14 

=ASG, T 

1 1 . ,F40/ l/POS/10 

is 

= A S G , T 

12 . , F40/1. /POS/1 0 

16 

= A S G , T 

13. » F 40 / 1 / PO S / 1 0 

17 

= ASG ,.T 

14. , F40/1 /POS/ 10 

18 

= A S G , T 

15. ,F40/ l/POS/10 

19 

= A S G , T 

1 6. , F40/1/PDS/1 0 

20 

= ASG ,T 

17 . , F40 / 1 /PO S / 1 0 

2L 

= A S G , T ' 

21., F40/1 /POS/ 10 

22 

=ASG,T 

22. ,F40/ l /POS/1 0 

25 

= A S G , T 

23., F40/1 /POS/10 

24 

= A S G , T 

24 . ,P40 / l / PCI S/ 10 

25 

= A S G » T 

25. , F40/ 1 / POS/1 0 

2,6 

= A S G , T 

26. , F 4 0 / 1 /PUS/ 10 

2? 

= A S G , T 

27. , F40 / 1 /PGS/1 0 

2% 

= ASG , T 

28 , U 

23 

= A S G , T 

30. , F40/1 /POS/1 0 

30 

= A S G , T 

31 . , F 4 0 / 1 /POS/10 

31 

= F L 1 , T L 0 DATA 

32 

RDM 1 40 

PHIL l PPUS 

33 

LA TER AL 

SLOSH FOR SPHERICAL 

34 

SPHERE 

RIGID 40 PERCENT 

35 

oooooooooo 

36 

I N I Tit. 


37 ■ 

G M X Y 7 


32 

X Y 7 FI )(_ 


33 

8.00 


40 

0.0026 2000000.0 -386.0 

41 

0.000228 400.00 0. 

42 

7 

12 

45 

AT 

l 5 

44 

1 

1 8.426723 

45 

1 

5 0.024364 

46 

oooooooooo 

47 

RP 

1 12 

42 

1 

1 1.00 

43 

1 

5 4.30 

30 

1 * 

9 7 . 3 

51 

OOOOOOOOOO 

52 

NP IS 

l 1 3 

55 

l 

15 4 4 


55 
5 i 
3? 


^JSSf 


' ANi; BLADDER 
# 386.0 IN / 


IN IT. 

S E C * # 2 . 


0.0 

n .06 


o.oo 


2.00 
6 . 40 
7.3 


0.0 


o 


1.164276 


, 73 
,20 
9 


- 0.336903 


3.30 
6. 9 
H.O 


ooooooooon 

ME LENT '] 
1 1 
oooooooooo 


12 
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58 

S3 

6c 

a 

6z 

63 

64 

65 

66 

67 

69 

70 

71 

7* 

73 

74 

73 

76 

77 

78 
73 
Sc 
Si 
$Z 
S3 
84 
Ss 
8<z 

87 

88 
S3 
5JO 

31 

32 

33 

34 

35 

36 

37 

38 
33 

IOO 

idi 

JOZ 

107, 

104 

105 
iC6 

107 

108 
103 
HO 
iti 

1LZ 

113 

114 

115 
116, 
ur 


0.0 

2.3 

4.1 

>.06275 

2.3 

4. 1 

0.3 

2.3 

4.1 

0.5 

2.3 

4. 1 

0.71 

2.3 

4.1 

1.3 

2.3 

4.1 

2.3 

4.1 

4.75 

2.95 

4.1 

4.75 

4.1 

4.75 

7 . 0 

4.7 5 

7.0 ■ 

9.0 

5.4 

7.0 

9.0 

7.0 

9.0 


8.00 




PEN TA 

T F T R A 

HEX A 

HFXA 

HFXA 

HEXA 

HEXA 

PENT A 

HEXA 

PENT A 

P F N1 A 

PENTA 

EjEXA 

PENT* 

PEW 1 A 

PENTA 

PE NT A 

PENTA 

HFXA 

PENTA 

PENTA 

PENTA 

F1NFL 


MODES 


IWMODE 

1 1 

0000000000 

NW 

6 

MU 

6 

SHIFT 

500.0 

MAX IT 

15 

MnPUNCH 


MFCHEQ 


X Y 2 

-2 PR UNI 40 

321 375 

XRFE 

1 3 

0000000000 

JVEC 

1 6 

1 

1 l 2 

0000000000 

JDUF 

-2 8 RUN 140 

EUL 

-28RUN .140 

MASS 

-2RR0N140 

MOOES 

-28RUN 1 40 

PLOT 


SURF AS 

1 13 

1 

1 33 26 

0000000000 

XY2 

-28RUN140 

EUL 

-28RUN 140 

FREO 

-2 8 RUN 140 

MODES 

-2 8 RUN 140 

jnoF 

-28RUM140 

S fART 


S TUP 


= END 


= XOT MM2 OOG 

=ADO,PL 

DATA 

= PMD » PLE 



6.68 

6.68 

5.13 

4.28 

4.75 

7.0 

7.0 

9.0 
10.2 
10.2 


0 3 


1.04 13 2 


7.57 


9.0 

10.2 

10.5 



ssssas? 


160 195 230 265 293 307 


314 
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lit 


= F I N 


\ 



